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ABSTRACT 


Tke primary intent of this work is to investigate 
computer aided compensator deSign using classical 
freguency respense technigues in conjunction with the 
technigues of modern mathematical prcgramming. A 
computer program to perform the automated design task 
is presented. This particular algorithm is based on 
the ccnstrained optimization technique introduced by 
M. J. ECX. 


In particular, the desired cpen lcop frequency 
response is specified for a number of discrete 
frequency feints over the frequency range of interest. 
Then tke ginimization routine is used to vary the 
compensatcr parameters in such a manner as to minimize 
a cost functional based on the difference between the 
actual and desired open loop frequency response of the 
compensated system. To illustrate the algcrithn 
several detailed examples using the program are 
Fresented. 
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I. INTRODUCTION 


Concurrent with the proliferation cf the digital 
ccmputer as a viable, indeed often indispensable, 
engineering tool, the development of state-space methods 
resulted in the application of the digital computer in the 
area of modern control theory over the past decade anda 
half. Initially this application was primarily cf an 
analysis nature, but as the theory of optinal control came 
into wider publicaticn and acceptance, the role of the 
digital ccmputer began to take on aspects of the éesign 
functicn. This situation is graphically illustrated by the 
myriad of industrial proprietary and university develcped 
ccmputer programs for computer aided control system design 
based ufecn state variable and optimal ccntrol thecries. 
However, as cointed out by Mitchell [11] comparatively 
little work has been done to extend and apply the 
ccmputaticnal fower, speed, and accuracy of modern digital 
ccmputers tc the proven classical theory and methcds of 
ccntrol system design. After this essential abandonment of 
the classical sethods (particularly by the academic and 
thecretical factions of the control engineering ccmmunity) 
in favor of the modern theory and techniques, recent 
puklicaticns indicate a realization of the need that both 
approaches to system design be used in a complimentary 
Cather than sutually exclusive manner. In particular the 
works of Ccffey, MacFarlane, Mitchell, Page, Rosentrock, 
Stear, [5,8,11,14] and others show a renewed interest in the 
tried and prcven classical thecry cf design coupled with the 
new technigues and approaches developed by proponents of 
mcdern ccntrcl theory, as applicable to the develcrpment cf 
usartle ccmputer aided control system design algorithms. 
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Ccntricl systen design based on frequency domain 
specifications, for linear time invariant parameter systems 
using classical methods, has traditionally been done 
essentially using trial and error techniques. That is, cnce 
the form and location of the compensator has been decided 
upcn, the design process reduces to cne of iterative 
analysis until the specifications are wet or it beccmes 
obvious tc tke designer that the compensatcr form chosen 
will not satisfy the specifications. When systems cf any 
complexity are invclved this method quite cften transforms 


intc an extremely tedious and time consuming process. 


The availability of a good analysis pregram coupled with 
a digital ccmputer operating in an interactive mode, to 
prcvide sinimal turn around time, can be used to somewhat 
Binimize the tedium of this type of design procedure. By 
varying the ccmpensator parameters the designer may then 
observe the effects upon system performance withovt having 
to expend ccnsiderable effort in repeated calculaticns. 
This design by iterative analysis may be ccntinued allowing 
the cosputer to perform the calculations until the designer 
is satisfied with the respcnse obtained. With the recent 
refinements in mathematical programming techniques, it is 
possible to relegate some of the design decisions directly 
to the ccaputer prcgram, thus even further minimizing the 
tiwe and efforts required of the design engineer. The 
ultimate coal, which is of course at present a long way down 
the road, wceuld be to input only the plant dynamics ard the 
desired specifications and have the computer perforgz the 
design process conpletely, without any interaction with the 
design engineer. 


In this thesis an investigation of the autcmaticn cf the 
design precedure is undertaken. In particular, giver the 
desired cpen loop response in the frequency domain a 
finimization routine is used in conjunction with an analysis 
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algcrithme tc systematically vary the compensator parameters 
tc approximate the desired response in acccrdance with a 
performance measure which indicates the deviation cf the 
actual open locp frequency response from the desired open 


iccp frequency response. 
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Ae GENERAL AUTOMATED DESIGN PHILOSOPHY 


The word design, as used in many engineering texts of 
various disciplines, runs the gamut of describing a process 
dealing essentially with scaling already existing 
configurations and components to one dealing with a kighly 
creative activity allowing a maximum amount of freedcm in 
satisfying a set of requirements or specifications. Ona 
Gay to cay Easis the majority of the design work 
acccmplisted Ly eéngineers lies somewhere Letween these two 
extremes, ccmbining elements of relatively reutine, straight 
fcrward calculations with those of creativity and 
imaginaticn. 


The cverall scope of the design process in general is 
indeed sc fcrmidable that it often defies complete 
descripticn in other than general terms, much less 
guantificaticn in precise mathematical language. This 
latter is unfortunately a necessity in achievine the 
autcmaticn cf this process within the limitations of the 
presently available hardware and software accessable Ly the 
engineer. Largely for this reason, when considerinc the 
automatic design problem it must by necessity be viewed, for 
the near future at least, in terms of a limited scope. That 
is, such items as creativity, experience, and intuition, 
which are undoubtedly elements of the total design process 
aS viewed ir its entirety, cannot be achieved by a cconzputer 
as the device exists at the present tine. 
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One cf tke key elements of any design prcecess is 
decision making. In the overall design process the engineer 
makes decisions based upon the information available to hin 
at the time cf the design undertaking. Quite often these 
decisions are based on a comparison of numerical values of 
various elements of the design configuration. It is ir this 
specific area of the overall design process that the rodern 
Gay computer can serve as a Significant tool for the 
encineer invwclved in design, particularly when a large 
number cf ccapariscns (expressable in numerical terms) are 
invclved. In other words, the computer can be used ina 
ccmkinaticn analysis and elementary decisicn paking role to 
aid the engineer in the overall design precedure. Ey means 
ef analysis cumerical values may be obtained, forming a 
basis frcm which comparisons can be drawn. The analysis 
portion is generally essential for exact design and forms 
the basis fer any optimization. The speed available in the 
digital ccmwputer allows many calculations and comfarisecns, 
over a wide range of possible solutions, to be performed 
which would Ee prohibitive if peformed by hand. 


The employment of the digital computer as an aid in the 
aggregate design process is perhaps as much an art ag _ the 
design process itself. While the use cf computer aided 
désign prcgrams can free the engineer from many burdensome 
and time ccnsuming calculations, the individual engineer 
must still previde the framework in which tc interpret the 
results in the hard light of physical significance. The 
engineer must also be capable and willing to define the 
Frcblem in a form with which the computer can work. More 
often than net, this latter requires a much sharper, clearer 
definiticn (and subsequently understanding) of the problen 
than if the same problem were teing designed long hand. 
Unlike the analysis problem, which generally possesses a 
specific answer, the design problem may not have a unigue 
Sscluticn which further complicates any computer aided design 
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algcrithz. Ancther aspect that should not be overlocked is 
the fact that when using a computer aided design program the 
engineer cften receives, if nothing else, an indicaticn of 
the soluticns which are not feasible. This in itself is 
often nct a trivial result especially in dealing with 
systems cf any complexity. 


Ancther important factor which should nct be overlooked 
is that cnce an answer is obtained using a computer aided 
design algcrithm scme type of sensitivity study will more 
than likely ke necessary. This is largely due to the fact 
that ccmpcnent values in phySical hardware cannot generally 
be contrclled to the precision with which a computer does 
its calculaticns. Therefore, while a particular algcrithno 
for use in automated design may generate what is indeed a 
theoretically flausible solution the implementation of such 
a soluticn gay not be possible if it is extremely sensitive 
te the parameter values involved. 


B. SOME ASFECTS CF COMPUTERIZING THE CLASSICAL DESIGN CF 
CCNTROL SYSTEMS 


Classically the purpose of the design problem in 
feedback contrcl systems has been to make a given plant 
satisfy a set cf performance specifications imposed upen the 
system. fLCesign of linear, time invariant ccmpensators for 
the ccntrcl cf linear systems has formed a major pertion of 
the design effcrts of conventional servomechanisms. In the 
process cf investigating the automation of this design 
precedure, certain key elements were found to be generally 
necessary in istpdementing a basic practical ccsputer 
algcriths for autorated design, regardless cf whether it is 
based on tine cr frequency domain approaches. These 
elements consisted of: 
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1. A ccmglete description of the plant 
cr process to be controlled. 


2. Specifications describing the 
response desired from the system. 


3. Some type of criterion functicn or 

Dem detdee ioe deeee os ace ae cae 

4. The form or structure of the 

it epsacor or ecntrol to be 
These items appear in varying degrees aS a common 
denominatcr amcng present day automated design algorithms 
{(1,5,11,13] regardless of the specific methods emfioyed or 
systems keing considered. That is to say, these elements 
appear in algorithms designed for time demain, trarsfornm 
dowain, stochastic frocesses, and nonlinear systems. 
Generally items 1. and 4. are employed in some fcrm of 
analysis algcrithm while items 2. and 3. form the basis for 
a Minimization (cr maximization) routine that prevides a 
foundaticn fcr some type of decision making withir the 
algorithn. This idea of using a performance measure upon 
which to kasé scme type of decision making routine is a 
direct carrycver from develcpments in the computerized 
design of oftimal control systems which relies heavily upon 
the idea cf minizizing some criterion function. Unlike the 
technigues cf cptimal control however, the cost functicnal 
as used in this ccntext is nct necessarily a function cf the 
states and tke applied control, rut rather a measure of the 
deviation of the actual response from some desired response 
which has been specified by the designer. This is perhaps a 
fine distincticn but one that must be recognized and kept in 
Bind. One of the nany consequences of this is that there is 
no ilicnger a strcng motivation to restrict the performance 
measure tc a guadratic form in order that the second sethod 
of Liapunov may te used to frove that a stable system will 
result. Rather, the form of the cost functicn in 
cCcmputerized design based upon classical ccntrol thecry is 


chosen tc have properties which aid in the convergence of 
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tke numerical mninimization routine used in the algcrithm or 
which has scre inherent physical significance in the design 


prcecess. 


C. A FCEMUIATICN OF THE AUTCMATIC DESIGN ERCBLEM IN THE 
FREQUENCY DCKAIN 


In this thesis the automatic design prcblem based on 
theory develcped for classical design methods using 
freguency dczain procedures is investigated. Several 
approaches FEased on frequency domain techniques have keen 
presented by various people to automate the control system 
design prccess [5,11,13,14,20]. These algorithms range fron 
methods cesigned primarily for use with multivariable 
systems to specialized interpretations of frequency domain 
data such as vector frequency domain technigues and inverse 
Nyguist setkcds. The suitability of these various methods 
have keen demcnstrated by their respective authors and 
others scrking in the area of computer aided control system 
design. While each represents a step forward in the line of 
progress tcward automating the design process using 
freguency dcmain techniques the Suitability of each 
algorithm is dictated to some extent by the particular type 
of design preblem being solved. Many cf these methods 
require gradient calculations, and those that do not rely 
primarily on simplex nonlinear programming algorithms. The 
former methcds restrict the ferm of the cost functicn to be 
Ririmized in the sense that it must be differentiakle and 
prcblems inherent with numerical differentiating algorithms 
are obvicusly involved. The latter methods invclving 
Sinplex algcrithms have difficulty in handling constrained 
Hinimizaticn problems, particularly when solutions exist 
hear or cn the constraint boundaries. 
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One cf tke advantages of the freguency response apfroach 
lies in the fact that it lends itself easily to the use of 
cpéen-loor frequency response data in designing closed loop 
systems. The use of frequency response methods provides a 
graphical pertrayal of the system response that lends itself 
well to quick interpretation Ey the designer. Also since 
approaches such as those developed by Bode were intended to 
give Ssitrle approximate methods for use in hand 
calculaticns, adaptation of these methods in a computer 
algorithms allcws the engineer te plan his design using such 
technigues and then employ the aid of the ccmputer without 
the requirezent cf significant problem restructuring in 
order to cast the problem in some particular format not 
generally familiar to the designer, but necessary for 
solution Ey a computer algorithn. With systems invclving 
ccmplex transfer functions such restructuring can itself 
become quite kurdensome thus defeating one cf the primary 
reascns for turning to computer aided design. Perhaps it 
shculd be peinted cut at this foint that even using a 
ccuputer algorithm to aid in the design process does not 
totally relieve the engineer of all hand calculations. Time 
must be taken to carefully map out a plan of attack cna 
particular freblem. In essence the designer must thrcughly 
understangd the problem and this is generally achieved by 
scgkeé preliginary rcugh sketches and calculations before the 


ccmputer ever sakes its appearance in the design process. 


Throughout this thesis, the vehicle used for presenting 
the autcmated design algorithm is the Singleé-input 
Single-output unity feedback system of the form shown in 
figure II-1. 
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General System, Block Diagram 


figure II-1 
The plant is considered to be linear, time invariant, and 
ccapletely kncwn. Thus the plant may ke descrited Ly a 
rational transfer function in the Laplace variable s of the 


form 
n 
ITs +z.) 
G(s) Zz i=0 1 
TT +ep.) 
2 
1=0 
where z and p. may be complex constants. The compensator 
1 a 


G (s) is alse considered to be a rational transfer furction 
c 


in s. 


The specifications are given in terms cf the system open 
locp freguercy response, and the desired respcnse is 
achieved by reshaping the system open loop frequency 
response thrcugh the use of the compensator in the _ feed 
forward fath. This basic procedure is well known tc most 
centrol engineers and is by no means a novel approach to the 
design roklem [3,.10,17,18], but through the years it has 
found extensive use in the practical design of ccentrol 
systems. 
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The methcd proposed and under investigation here is a 
relatively straight forward adaptation of this technigue 
with a ccmuputer algorithm harnessed to do the ccmpersator 
parameter variaticn in order to systematically reshare the 
system ofen locp frequency response with a pinimum amount of 
calculaticn and iteration required on the part cf the 
designer. The desired open loop response is specified by 
means of a sagnitude and phase profile cver the freguency 
rarge of interest. This is nothing more than a set of 
desired magnitude and corresponding phase values at discrete 
frequency pfeoints selected by the designer. In addition a 
specific form or structure is postulated for the compensator 


transfer function G (s). At this point cf the autcmated 
c 


design Erocess the algorithm relies heavily cn the 


experience and irtuition of the individual designer in 
checsing a structure or form of compensatcr that will have 
scmeé reascnakle prcbability of success in reshafinc the 
frequency respense curves in order to satisfy the design 
specificaticn. Since no firm analytical kasis exists for 
the deterrination of a Specific structure for the 
ccmpensatcr, the computer algorithm must be given this 
informaticn or an iterative method based on some type of 
criterion must be used in order to estaklish the form. This 
latter swethcd was felt to ke much too time consuming and 
restrictive, because of the wide range of variations in 
Fessible structures that might have to ke considered, and 
essentially reyond the scope of this investigation. Thus 
the choice of the compensator form is relegated to the 
program user as a design decision that must be undertaken 
Prior tc using the algorithm to determine specific values 
for the ccefficients of the compensator transfer function 
that will yield the desired response. This frocedure should 
net present a major stumbling block for the designer 
prcvided that he possesss a thorough understanding cf the 
prceklem tc te solved. In fact, by inputing the specific fcrn 
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of the ccmpensator the designer may investigate compersator 
configuraticrs which might otherwise be considered 
unfeasible fer the solution of the problem if some type of 
iterative algcritha were used in determining the structure. 
With these items ferming the basic inputs tc the algcrithm a 
ncnlinear perfcrmance measure coupled with a minimization 
reutine is tken used to minimize the difference between the 
desired and actual frequency resfronse of the ofen loorp 
transfer function given by, 
T(s) = Cate ee 


The perfcrmance measure (or cost function) selected is cf a 
normalized form frimarily as a matter of convenience in 
interpreting the algorithm's results regardless of the 
particular jyreblem being considered. That 1s , the cost 


function kas a general form given by 


A (jx) 
D (jy) 
where D(jw) represents the desired open loop frequency 


J(jw) = £(1.0 - 





response cf the entire system, A(jw) is the actual ofen loop 
response, ard w represents the set of discrete freguency 
values at which the cost function is to be evaluated. AS 
mentioned previously this general form cf normalized cost 
function was selected primarily as a matter of convenience 
in that a farticular cost function can be constructed in 
such a manner that the optimal or smallest value cf J (jw) 


under ideal scluticn cf the problem will be zero. 


By inguting the desired response in the form of the 
magnitude and phase prefile described above, ncrmal 
frequency dcnzain specifications such as gain margin, fhase 
Margin, and tandwidth (here defined as the open loorg gain 
crossover frequency) may be supplied to the algorithm in one 
type of fcrmat. With these specifications forming part of 


the gairc ard phase profiles they are automatically 
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inccrporated into the overall cost functicn each tise the 
algorithm e€veluates the frequency response cf the ofen loop 
System cver the range of frequencies specified ky the 
designer, thus also eliminating the requirement of 
specialized ccmputation within the algorithm to check for 
satisfacticn cf these particular specifications as the 


compensatcr farameters are varied. 
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A. GENEBAL 


A variety of functional minimization techniques have 
ccue into use in recent years are discussed in ref. [8]. 
The minimization algorithm used in this work was developed 
originally ky M. J. Eox, [4], and is capable of finding the 
hinimum cf a general nonlinear cost functicn, composed of 
several variatles, within a ccnstrained regicn. The ccmplex 
method of Eox has reen implemented in program form as_ fart 
of the stancard sukroutine library at the Naval Postgraduate 
School W&. EF. Church Computer Center. The subrcutine 
(BCXPLX), was originally programmed by R. R. Hilleary cf the 
staff of tke Naval Postgraduate School Computer Center and 
has subsequently been modified slightly Ey the author in 
order to render it more suitable for use in the scheme of 


the compensatcr optimization progran. 


In additicn to its availability as a standard 
subroutine, several other features of this farticular 
finimizaticn algorithm made it an attractive choice fcr use 
in minimizing a cost function based on frequency domain 
specifications. While admittedly more efficient and 
sophisticated algorithms for function minimization are 
available, tke flexibility and ease of programming of the 
cctplex eetkcd of Bcx are among considerations that should 
not be overlccked. This method does net require that the 
derivatives cf the cost function or implicit constraints be 
calculated, thus avoiding the sundry pretlems typically 
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asscciated tith numerical methods for the calculaticn of 
derivatives, and eliminating the restriction that the 
derivative cf the cost function exist over the range of 
interest. In implementing the minimization method the cost 
function need not be expressed explicitly in terms cf the 
variables tc be manipulated by the program in the 
Birimizaticn process. Unlike many of the cther ninimization 
prcgrams avéilable, BOXPLX, also provides a restart 
capability. That is, once a possible minimum value cf the 
cost functicnal is reached the program automatically 
restarts frem rardomly generated initial values cf the 
variables. f&hile this method by no means guarantees tkat a 
glcbhal ewinimum will be achieved for every conceivable 
prcktlem the restart capability does at least sorewhat 
attempt tc avoid the problem cf the minimization terminating 
on a local winima within the search area. This methed is 
also tc scme extent scale independent in that the size of 
the initial cecmetric figure used in the search for the 
Binimum value cf the objective function is rcughly scaled to 
tke crder cf the problem variables, through the use of the 
difference ketween the lower and upper bounds cn the 
variables. 


As is the case with other sinimization algorithms, this 
Particular cre is not void of inherent disadvantages and 
difficulties. Present in any minimizaticn algorithm of 
course is the recurring problem that no guarantee can be 
made of achieving a global minimum fer all classes of 
problems. The ccmplex method cannot handle eguality 
constraints without modification cf the algcrithm. Also as 
the number of variables increases the method rapidly becomes 
inefficiert. Finally the unconstrained probkleng is 
apparently mcre efficiently handled via gradient technigues. 


B. DESCKIFIICN OF THE ALGORITHM 
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The giniwization method developed by Ecx has its kEasis 
in the seguential simplex techniques of linear prograagnming, 
but is specifically designed to obtain solutions to freblems 
cf ccnstrained optimization and to avoid such prokless as 
the inakility of ccnstructing the first siuaplex. Thus the 
resulting scluticn of the minimization problem is a set of 
variables wkich may lie at the extreme edges or boundaries 
cf their permissible ranges. The proklems generally 
associated with ccnstraints on the variables is essentially 
attacked Ly tsing a flexible geometric figure having n+l, or 
more, vertices and capable of expanding or contracting, in 
any or all directions. Therefore, unlike many of the 
Sinplex alccrithms no provision is made to maintain a 
regular cecnetry between the various ” points in an 
n-dimensicnal space, but rather a changing structure (dukbed 
a ccmplex), capable of flattening, rounding corners, and 
eventually ccllapsing on the point representing the ninimun 
of the oktective function, is employed rendering the sethod 
more flexible in handling constrained minimization froktlems. 
The numker cf independent variables in the cbhjective 
function and ccnstraint equations dictate either directly or 
implicitly tke dimensionality cf the space tc be searched. 
The vertices saking up the complex are continually rejected 
and generatec as a search for the minimum cf the okjective 
functicn is carried out with each new vertex generated being 
reguired to satisfy all the imposed constraints. 


The ccnmplex method searches for the minimum of an 
objective function J(x), where the vector x represents a set 
of n variatles, within a region of space Ecunded Ly ufrer 
and lower linits on the variables (explicit constraints) 
given by, 


ip SU, «61=1,2,3,;-.-,0 
1 al 1 
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and in addition subject to the restrictions imposed by the 
ccnstraint furctions (implicit constraints) cf the forn, 
1, @ 2 0 J=1,2,37000 70 


As with tke simplex methods this is an iterative algorithn. 
It uses a set of k 2 n+1 points which form the vertices of 
the complex and simultaneously satisfy all the isgposed 
ccnstraints. Initially only one pcint or vertex oust 
satisfy all the comstraints and serves as a means of 
generatirg the other k-1 vertices through the use of 


pseudorandce rumbers r,., uniformly distrirkuted cver the 
4 


clesed interval from zero to one. The additional k-1 


vertices cf the first complex are established by the 
fclicwing fcraula: 

x =L +r e(U. - 1.) 

1 ah i a a 
As can re seen frem the formula, the vertices generated in 
this fashion will always satisfy the explicit constraints, 
however, €ach cne pust be checked to insure compliance with 
the implicit ccnstraints. [f) an aoplicrt constraint 
viclation is indeed found to cccur at this stage then the 
trial vertex is simply moved, repeatedly if necessary, 
halfway in tcward the centroid of the cther already accepted 
vertices, until ultimately a feasible eint is fcund. 
Continued application of this precedure will result in the 
gceneraticn cf the k-1 additicnal points required fcr the 
initial ccmplex within the feasible region defined by the 
ccnstraints. 


Cnce the initial complex has been determined the 
iteraticn preceeds ty evaluating the objective functior J (x) 
at €ach vertex. The vertex at which the oktjective function 
has the largest value is then designated as the current 
worst vertex and this point is reflected througk the 


centroid of the remaining vertices establishing a new 
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W 
complex. If this worst vertex is designated as x then the 
overreflecticn is accomplished according to the formula 


N CR W 
xX = (1.0 + alpha) ex - alphaex 


CR 
where alpha 2 1.0 is the overreflection ccefficient, x 


Lepresents the centroid of the remaining vertices,and oa 1S 
the new vertex. If the overreflection of the current werst 
vertex resuits ina violaticn cf the constraints or if the 
new point still has the worst value in the set of vertices 
then the feint is mcved halfway toward the centroid used in 
the overreflection process. This retraction is repeated (if 
necessary) until the overreflection coefficient is reduced 
to some stall value, delta, or a Suitable new vertex is 
established. If the overreflection coefficient is reduced 
to the value delta and the objective functicn value at this 
new foint is still the worst in the set of vertices then the 
prcejected vertex is replaced by its original value and the 
second wcrst vertex will be overreflected instead. This 
prccess keeps the complex moving toward the winimumn of the 
ckjective function provided the figure has ‘not ccllapsed 
mato its centroid. Once assured that the new ccmplex 
Satisfies all the constraints and results in an improvement 
of the okjective function, vertices of the complex are 
rejected and generated in a systematic manner aized at 
Binipization cf the objective function. Thus essentially 
the complex moves over the feasible region eventually 
straddlin¢c and collapsing upon the point that renders a 
Birimum value for the objective functicn. 


This frecess ccntinues and is eventually terminated when 


the complex shrinks to a predetermined acceptaktle small 
Value epsilcr given by: 
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k 
1 Cc v 2)0.5 
‘ - Dt sie) - 3) ] < epsilcn 
1=1 


c 
Where J(x ) represents the value of the okjective function 
V 
at the centrcid ard J(x ) is the value at the varicus 


vertices. 


As pcinted out by Box [4], the exact values cf the 
overreflecticn coefficient, alpha, and the exact numter of 
vertices, k, do not appear to be critical to the functioning 
of the algcrithm prcevided that they are greater than unity 
and n¢t1 respectively. The main concern here is that a value 
cof alpha larger than unity keeps the ccmplex from shrinking 
prematurely as it traverses the feasible region in search of 
the pinizum ckjective function value. Values of alpha = 1.3 
and k = zn as given in ref. 4 are used here since the exact 
values chosé€r are apparently not critical and good results 
have been reperted, by the originator of the algorithna, 
using these farticular values. For a detailed example of 
the iteraticn process and a more sophisticated explaration 
of the algcrithm the interested reader is directed to refs. 
2, 4, and 6. 
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A. DESCEIETICN OF ALGORITHS IMELEMENTATICN 


The frequency domain algorithm described in the latter 
part of Chapter II has been irplemented in program forn 
using FCRTEAN IV ccding and ccnfigured for use on the IBM 
360-67 ccrputer system available at the Naval Postgraduate 
Schcol Ccmputer Center. With the exception of the graphing 
sukroutine tke remainder of the program is self ccntained 
and may be used on any computer system having a standard 
FORTRAN IV compiler. Without the overhead routines 
reguired of the graphing subroutine, the basic prcgram 
reguires apfroximately 180k bytes of ccre memory, and 
including tke core requirements of the flotting routines 
roughly 210k ktytes of memory are reguired ky the prcgranm. 
While admitedly less core would be necessary if such 
techniques as array overlaying were employed when 
prcgramming the algorithm, it was felt that for sinuplicity 
in following the flow of the program by future users these 
technigues were best left unused. 


The ccmputer aided linear time invariant compensator 
optimizaticn pregram (CALICO) is essentially composed of 
fcur main farts with numerous additional small sukroutines 
functioning in supfortive roles. Since often a progras user 
May wish tc ferform certain modifications tc programs of 
this nature in order that they may better serve his 
particular purpose, Appendix A provides a descripticn cf the 
SUffortive sukroutines used within the frogram to facilitate 
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acccmplishment of any program mcdificaticns that may be 
desired. A block diagram of the major components cf the 
overall program are shown in figure IV-1. The main program 
serves frimarily tc input data describing the system under 
consideraticr and to output the values of the roots and the 
gain of the compensator that will best achieve the desired 
response in accordance with the cost function selected. 
During program execution, the main frogram also serves to 
@Girect the seguence of execution of the other major sections 
of the CALICC algcrithno. 


Since the plant and compensator are assumed to be linear 
time-invariant, they may be completely described by transfer 
functions in the independent Laplace variable s. Once the 
plant and assumed f¢rm of the series compensation are input 
alcng witk the profiles forming the desired phase and 
magnitude respcnse program control Froceeds tc _ the 
Elnimization subroutine, BOXPLX. From the point of view of 
computer aided ccntrol system design, asscciated witk this 
sukroutine are two key function sSub-programs, FE(X) and 
KE (X). The function Sub~program FE(X) performs the 
calculaticn cf the system open loop frequency response, at 
the discrete frequencies provided in the input data, and 
calculates tke value of the cost function FEased upon the 
difference tetween the actual and desired frequency response 
at the discrete frequency values. Six separate cost 
functions, frem which the program user may select cne, are 
incerperated into this function sub-programp. The value of 
the cost function is used by the ninimization routine to 
vary the ccmpensatcr parameters in order to minimize this 
quantity. The function sub-program KE(X) checks fcr any 
Violaticn cf the implicit constraints, which in this 
particular case are concerned mainly with the presence of 
right half plane rcots. The iaplicit constraints used, also 
impose limits on the real and imaginary parts cf the 
nhuperator anc denoginator polynomials evaluated at any 
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MAIN PROGRAM 


Output optimized results; call simulation, 
Sitebiltty, and miotting nemtines - 


SUBROUTINE BOXPLX 
Systematically varies the compensator 
parameters to minimize the cost function J. If 
the number of trials is exceedid the best 
solution achieved to that point is selected as 
Mie sOluUtton. ££ tne cost function remains 
unchanged for a soecific numoer of trials the 


minimum is considered as having been achieved. 


PUIG TLEON FE PUNCT Ee mks 
Calculates the fre- Checks fOr sigugnt 

suency response of nalf olane roots 
the open loop system and violation of 
at the discrete val- Ocher implreve 
ues svecified in the constraints. 
input and forms the 
difference between 
Geslwed mand ac cual 
response making up 
icwee COsemeunction. 











SUBROUTINE SIMULA; Calculates the open loop | 
frecuency response using ontimized parameters. r< 


SUBROUTINE GRAPH; Plots the desired graphs 
Svcecified bv the user. 


Program Block Diagram 


Figure IV-l 
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particular frequency. These latter limits are a direct 
result of the finite computer word size and the particular 
methcd tsed to represent numerical values withir the 
machine's mencry. Bcth of these factors have an effect in 
deterpining the saximum and minimum size numerical value 
that may ke represented within any fparticular comfuter. 
While the ccefficient values determined by the minimization 
routine will nct exceed these waxinum or pinimum values, the 
real and isaginary parts of the numerator and denominator 
EPelynomials evaluated at a particular frequency are squared 
during the frocess cf determining the tsagnitude cf the 
frequency respense. Without the implicit ccnstraints these 
sguared values can exceed the size of the numerical 
capability of the sachine, particularly in the case cf high 
crder systems, causing premature termination of the 
algcrithm ané hence failure in obtaining a_ solution. Cnce 
tke miniwizaticn subroutine determines the ccefficients that 
result in a ninimum value of the cost function cr the 
Ginimization process is terminated for other reascns, 
prcgram ccntrcl is transfered back to the main fpregran. 
Here the reots and gain of the compensator which resulted in 
a GBinimum cest function are calculated and cutput. The cpen 
locp fregquercy response of the sytem is then calctlated 
using a mcre densely distributed number of frequency values 
than would normally be selected by the designer fcr use in 
the minimizaticn process. This is done largely to show the 
existence cf any resonant feaks which might otherwise be 
straddled and not apparent if only the discrete frequencies 
selected Fry the designer were used in the plotting routines. 
The stability cf the entire system is then checked by adding 
the open lccp numerator and denominatcr pelynomials and 
applying the Rcuth stability critericn tc the resulting 
FOlynomial which is the system characteristic equation. 
Firally the results at the discrete frequencies specified by 
the designer and the "continuous" curve generated Ly the 


SIKFULA sukrcttine are output cn Bode and gain versus _ fhase 


a2 





Picts. In addition, plots of the difference between the 
desired and actual results at the discrete frequencies 
specified by the designer are output in crder to aid the 
designer in sgaking any modifications to the form of the 
cc&pensator that may be necessary in order tc better satisfy 


the system specifications. 


As stated in the previous paragraph six separate cost 
functions have been incorrforated into the furction 
Suk-program FE(X). The first three of these are cf the 
€Xrror squared type and the latter three are cf the aksclute 
error type. These six cost functicns represent various 
ccmbinaticns cf one minus the actual response divided Ey the 
desired respense at the discrete frequencies under 
consideraticr. These pre-defined cost functions are 
prcevided as a matter of convenience in using the prcegram and 
certainly de nct prevent the user from defining his cwr cost 
function previded this is placed within the FE(X) function 
suk-program. The first cost function combines the deviation 
of koth tke desired magnitude and phase resfonses frem the 
actual magnitude and phase responses of the ofen Ilcop 
system. That is the first cost function is cf the fors 


I) “] («= eae : (1.0 P i 
DMAG (a) DPHASE (o) 
where AAG and APHASE are the actual magnitude and phase 
response cf the cpen loop system and DMAG and DPHASE are the 
desired sagritude and phase respectively. These quantities 
are evaluated in the cost function and summed cver the 
discrete frequency values in the specified range. The 
second ccst functicn which the user may chocse, uses solely 
the deviaticn of the actual magnitude response desired 
Magnitude response evaluated at the discrete frequency 
Values. Thus this cost function is of the general fers 
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( AMAG (w) ) 
J (w) = 1.0 - ——_._____ 
Sy DMAG (w) 


If a winisum phase system is being considered the magritude 
and phase fylcts as a function of frequency are uniquely 
related ty the Bode theorems {3]. That is, determining the 
Magnitude ratio cf a system over a given frequency range 
will alsc resvlt in a unigue phase relaticnship for the 
system. Thes if the systen design being considered is of 
the mininum phase tyge the use of the type two cast function 
described skculd result in both the magnitude and phase 
specificaticrs being satisfied. The third pre-defined cost 
function uses tke deviation of the actual open lcop fhase 
response fcrs the desired phase response in a functicn of 
the form | 


J (w) >> (1-0 ;, eee 7 
DPHASE (w) 
While for swinisum fhase systems the magnitude and phase 
responses are uniquely related, it should te noted that the 
fixed lcss cr gain of the system cannot be determined from 
the phase characteristics alcne. Consequently the use of 
the type three cost function may result in a constant gain 
errer fcr tte magnitude response of the system. This error 
however is readily observable from the output produced by 
the program and correction of this may be accomplisked by 
adjustment of the compensator gain value accordingly. The 
cther three pre-defined cost functions are of the same 
general fcrgz as these just presented with the exception that 
the squared differences are replaced by the absolute value 
of the differences. The general rationale in selecting this 
latter type cf cost function was that on occasicn the 
squared type cost function will have a tendency to flatten 
near the gininmum. With the absolute value ccst functicn the 


user may avcid the excessive time required in findinc the 
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Binimum valtes wten the case arises where the squared cost 


furcticn Leccmes flat. 


In utilizing this design program the user Should alsc be 
aware that the need fer frequency scaling cf the problem may 
arise. This is particularly true when high order systems 
involving high frequencies are being considered. The need 
for freguency scaling is primarily a result of the finite 
word size available in the computer system. It must Le kept 
in mind thet the magnitude response of the system is 
computed as the ratio of two ccmplex numbers and while the 
resulting macrnitude may be relatively small the value cf the 
Bagnitudes of the individual complex numbers forming this 
ratio may te large. The situaticn is further complicated by 
the fact that in ccmputing the magnitudes of the ccmplex 
nhumkers invclved the real and imaginary parts must be 
sguared, thus effectively reducing the maximum power of the 
numbers tshich can be allowed by a factor of two. If the 
magnitudes of the resulting numkers kecome too larce an 
errcr messace will ke printed recommending that the problen 
be freguercy scaled. 


Be DATA INECT ANT OUTPUT 


In crder to sisplify use cf the program, an attemft has 
been made tc maintain the input information necessary for 
executicn cf the algorithm to a minimum, and at the same 
time preserve a certain amount of flexibility withir the 
Frogram which may be controlled by the user. In so doing 
scme of the frogram control variables for various cftions 
have been set irternal to the program and may not be 
Manipulated ky the user in the input data deck. These 
internally set variables do not detract from the general 


functionirg cf the program, but mention of this is made to 
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bring it tc the user's attention in the event that 
specialized rfrceblens might more effectively be investigated 
ky resettirg scme of these normally preset prcgram 
variables. Many of the preset items are concerned with 
centrcl cf the tinimizaticn routine. Such items as the 
upper and lcwer bounds of the search area, the integer 
prcgramming cption, and the seed for the randcm number 
generator are set to specific values within the main frogran 
just pricr tc calling the minimization subroutine. For 
example fcr the ginimum phase case, the Ilcwer and upper 
bcunds cn the search areas of the compensator numeratcr and 
dercminatcr pelynomial coefficients are set at 0.C and 


6 
1.€x10 respectively. For the non-minigum phase case the 
beunds on the seé€arch areas cf the numeratcr pclynonial 
6 6 
coefficients cnly are extended to -1.6x10 and 1.6x1C and 


the implicit ccnstraijnts are adjusted to allcw fright half 
Flane rccts for the numerator polynomial. The integer 
prcgramming cfpticn for BOXPLX is also preset:'to a value of 
z€xro indicating that values for the pcliynomial coefficients 
need not ke integers. In general, selection of the irteger 
prcgragming cfeticn results in considerably more execution 
time fcr the program and for the purpose cf varyinc the 
compensatcr sarameters the selection of this option precvided 
ncthing tc erhance the soluticn of problems. 


The input data reguired for program executicr and 
control may te divided into general areas as follows: 1.) a 
title card; 2.) a series of cards giving the gain and 
ccefficierts cf the plant transfer function; 3.) a series 
of cards giving the assumed gain and coefficient values of 
the compensatcr transfer function; 4.) a control card 
describinc the flots desired, the type of cost functicn to 
Fe used, whether the system contains a zero crder hold, and 


tke number and range of the discrete frequency points being 
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considered; £.) cards containing the discrete frequency 
values, the desired magnitude, and the desired phase; 6.) 
finally a card describing the number of iterative trials to 
be allowed in the minimization and whether any cutgfut is 
desired frem the minimization routine in order to keep track 
of the sganner in which the coefficients are being varied. 
The specific input formats of the above necessary data are 
described ir detail in the following faragraphs. These 
paragraphs are numkered so as to indicate into which cf the 
six general areas a particular data card belongs. This 
descripticn ccupled with the program listing given in 
Appendix EF shculd frovide sufficient informaticn for setting 
up the data deck necessary for execution cf the fregran. 
Figure IV¥7r2 illustrates the data deck used in executing the 
first example precblem presented in section C of this 
chapter. In this figure each line of data represents an 
incividual data card. The column numbers for the data cards 
are alsc shcwn abcve and below the data deck secticn cf the 
figure. 


eS —- << 


4d. Title card. 


This must be the first card of the data deck and 
Frovides a s#eans of problem identification if several 
sé€parate prektlems are rcun under one jok name. Columns 1 
through 4E may be used to input whatever alphax~<numeric 
characters the user desires for identifying any particular 
Froklem. This title will also he printed cn the grafghical 
output reguested by the user. 


Iho 


A, Flart cain 


The seccrd card contains the gain of the plant in 
cclumns 1 through 10. If the gain is not explicitly present 
in the plant transfer function then the value input cn this 
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card shculd te unity. The format used tc read the plant 
gain is ar £10.00 type format. Thus, the gain may be input 
in either floating point or scientific notaticn type 


formats. 


This card specifies whether the numeratcr of the plant 
transfer functicn is to be read in factered cr pclyncnial 
form and what the crder of the numerator is. A P or an F in 
cclumn 1 indicates whether the plant numerator data is tc be 
read in pclyrcmial or factored form respectively. Columns 2 
and 3 cortain the crder of the numerator. Thus an input of 
the form FC3 starting in column 1 indicates that the 
EClyncmial is of third order and is to be read in factored 
fora. 


2C. Plart transfer functior numerator 


The fclicwing card(s) contain either the pclyromial 
coefficients cr the factors of the plant transfer function 
nugwerator, depending on which form has teen chcsen for 
aoput . If the pclynomial form of input has been selected 
then the gcelyncmial coefficients are input in ascending 
Powers of ¢ using an 8£E10.0 format. If more than eight 
coefficierts are needed they are simply continued on 
successive cards. The coefficient of the highest order tern 
is assumed tc be nermalized to unity and is not read in 
explicity. Thus the number of coefficients actually read 
ccrresponds tc the order of the numerator folynomial. rf 
the factcred forms of input is chosen then each factor is 
read in cn a separate card with columns 1 through 10 being 
used for the real part of the factor and cclumns 11 thrcugh 
20 for the imaginary part. Again the format used is £10.0, 
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thus allcwirtg the factors to be input in either flcating 
Feint or scientific notation forn. Since complex factors 
will occur in ccnjugate pairs the values cf these factors 
need only te input cnce. Thus the number cf cards used will 
vary deperding upon how many purely real and how many 


cormplex factcrs are involved. 


Sooo on = a a a _— Se ee 


Again tke fcrm of input is specified as in card number 
22 along with the order of the plant transfer furction 
dercminatcr. 


Ino 
lr 


lant transfer function derominatcr. 


ee eae a a SS a = Se a ee ee Se 


The fcrmat used to input the denominator is identical to 
that used fcr input of the numerator, contingent upfen the 
desired fcrg specified for the input. That is to say, it is 
possitle tec input the numerator in polyncnial form and the 
dencminatcr in factored form, or vice versa, since each tine 
a folyncsial is input as data the input form is uniquely 
specified. 





Folicwing the input of the information about the plant 
the descrifticn cf the form cf the compensator must be 
Supplied. This is begun by inputing an assumed comrersator 
gain using tke same format as was used to specify the plant 
gain. Cue to the manner in which the compensator transfer 
function nugerator coefficients are handled withir the 
Frcgram, if the assumed compensator gain is set to zerc, for 
lack cf any retter value, whis will result in the starting 
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guess of all the numerator coefficients being zero. fMhis 
will not cause the program to terminate prematurely, but the 
initial suess of zero for all the coefficients may not be 
what is desired by tke user. If doubt exists as to what 
value gain shculd te assumed, a value of 1 is suggested as a 


convinient arkitrary starting point. 


gut form and crder of the assumed ccmpensatcr 





The sethcd for specifying the form of the input fcr the 
ccmpensatcr initial numerator parameters is identical as 
that used in specifying the plant numerator. A P indicates 


eclynomial fcrm and an F indicates factored forn. 


In crder to maintain continuity in the input fcrmat the 
scheme and fcrmat used to input the compensator numerator 
parameters is identical to that used in specifying the flant 
transfer functicr numerator. 





As in previous cases either polyncemial cr factcrec¢ form 
is specified along with the crder of the compersator 
dercminatcr. 








Once again the same format as used to input the flant 
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éata is used in order to maintain the ccntinuity cf the 


foeut forsatirg. 


4a, Pregqram ccntroi card 


—— 


fhe Erecram control card is used tc input various 
quantities ccntrclling both the range and number of discrete 
peints used in the minimization of the ccst functior, the 
output desired frcw the program, and also whether a zero 
order hcld is present in the error channel of the system. 
As many as twelve separate quantities may be input on this 
sirgle ccntrel card in order to direct executicn cf the 
program. Cclumns 1 through 10 specify the minimum frequency 
(WMIN), in radians, of the total range Lteing consicered. 
Cclumns 11 through 20 specify the maximum value (WMAX), in 
radians, cf the range of freguencies being considered. Poth 
of these values are read using a 2E10.0 format. Fcllicwing 
the waxisun and minimum frequency values, up to ten input 
quantities are read using a 1013 format. Columns 21 through 
23 specify the numker of discrete frequency foints (NOMEG) 
over the given range at which the open loop frequency 
response is to be determined and the cost furction 
evaluated. The next input quantity (KNOW), in columns 24 
through 26 specifies whether the discrete frequencies of 
interest are to be read from data cards or automatically 
computed. Ncrmally this quantity is set to unity indicating 
that the discrete frequency points of interest will be read 
frem data cards. A value of zero indicates that the 
discrete frequency values are to be incremented linearly by 
the program starting at the minimum frequency value. A 
value of two indicates that the discrete frequency feints 
will be incremented logarithmically begining witk the 
tinimum fregvency value. Columns 27 thrceugh 29 indicate 
whether a Ecde flot of the results is desired. If this 
variable (NECCDE) is set equal to one no Bode flot will he 
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output. A value of zero for NBODE will result ir EFoth 
Bagnitude ancé phase flots. Columns 30 threugh 32 ccntrol 
the output cf the Nyquist plot. If NYQST is set to unity 
then the Nyguist plot will ke suppressed. A value cf zero 
for NYQST will result in a Nyquist plot output, however this 
is plotted crly over the range of frequencies between WHIN 
and WMAX. Cclumns 33 through 35 indicate the presence of a 
zero crder hcld (1Z0H) in the error channel of the unity 
feedback system. With IZOH set to zero the system is 
treated as keing of a continuous nature. With IZOE set 
equal tc cne a zero order hold is considered to be present 
in the error channel and the open loop magnitude and fhase 
computaticns are modified acccrdingly. Cclumns 36 through 
3& (NICHCL) specify whether a magnitude versus phase plot of 
the results are desired. As with the other flot opticns, a 
value of zere will result in a plot being preduced while a 
value of cre will suppress the magnitude versus phase flot 
as part cf tke output. Columns 39 through 41 contain the 
variable indicating whether non-minimum phase soluticns are 
tc be considered. With this variable (NMINFS) set equal to 
zerc only sinimum phase solutions will be allowed. A value 
of cne indicates that a non~minimum phase sclution will be 
allowed provided that it is the one which rinimizes the cost 
functicn. Cclumns 42 through 44 indicate if the desired 
magnitude prcefile values tc be input are in decibels. If 
this variable (IDB) is zero the magnitude values to te read 
are diszensicnless gain values, while a value cf one 
indicates that the magnitude values to be input will te in 
decibels. cCclumns 4S through 47 contain the variable IPLOT. 
This guantity is used to choose between printer flots or 
Calccmp flets for the graphical output. A value cf zero 
will result in printer plots and a value of unity will 
previde cCalccmp cutput of the graphs desired. mt ois 
reccmmended that the printer plot opticn be used in 
initially sclving a particular problem since this type of 
output recuires less turn around time than the Calcomp 
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output. Finally columns 48 through 50 specify the tyre of 
cest function (ICOST) the user wishes to chcose amoncg the 
six preprogrammed functions available in subroutine FE(X). 
The value specified should be between 1 and 6 correspcnding 
te the farticular cost function desired. If this value is 


zere the procram will default to the type cne cost function. 


Sampling period. 


This card must be inserted only if a zero order hcld is 
specified in the error signal channel. The sampling feriod 
T, in seccnds, is specified in cclumns 1 through 10. If the 
system kéing ccnsidered is of the continuous type then this 


card must ke critted. 


Next fcllcws a set of data cards specifying the values 
in radians fer second of the frequencies at which the 
desired gain and phase values of the open loop frequency 
response are specified and at which the cost functicn is to 
be calculated. These values are input in an 8£10.0 format 
using as many cards as necessary to speciry all the 
frequency values at which the program is to execute the 
alccrithg. 


This set cf data cards contains the desired ofen loop 
gain respcense at the discrete frequency values previously 
specified. MThe same 8E10.0 format is used here as was the 
case in reading the discrete frequency values. The same 


nugker cf cards shculd be used as was used in specifyirg the 
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discrete frequency values. 


Here fcllcws a series of cards specifyirg the desired 
open locg fhase response, in degrees, at the discrete 


freguency values being considered. 


Here toth the number of trials to te allowed in the 
Zinimizaticn frocess and the print interval for diagrestic 
purpeses are specified. These quantities are input using a 
215 format. If the number of trials is zero then a default 
value of 2C€0C is assumed. If the print interval is 
specified as zero then no output from the minimization 
routine will result. Caution should te exercised in 
Sgecifyinc the print interval in that if this interval is 
made too snall an excessive amount of printed output may 
result. The frogram user is referred to the comment section 
of subroutine BOXPLX in Appendix A for guidelines in 
selection of the print interval. 
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Figure IV-2 


Cc. SERIES CCMPENSATOR DESIGN EXAMPLE FPROELEMS 


In this section a number of example problems are 
presented and discussed in order to illustrate the use cf 
the computer aided compensator design algcrithm. In the 
hepe cf presenting a clear, concise illustration of the 
algoriths and its use the example problems begin with 
relatively simple straight forward text bock examples and 
gradually regress tc problems of a slightly more difficult 
nature. Alsc included are problems that are intended to 
illustrate scme limitations of the program in that these 
prctlems were not completely solved by the algorithm cr they 
may have xreguired a certain amcurt of manipulation cf the 
input data in crder to obtain a_soluticn. These latter 
types of recklems have been included since it is impcertant 
that the ¢esigner using any computer aided design algcrithm 
understand the algorithm's limitations and peculiarities as 
well as its ixtended capabilities. It is nct intended to 
gislead the reader by implying that the example frecblems 
presented here provide a comprehensive presentation of all 
Fossible Jlisitaticns of this automated design prcgran. 
Rather ttey represent those limitations which were 
enccuntered ir the process of verifying the program cver a 
necessarily linited period of time. With this in mind the 
first example precblem which deals with a relatively simple 


lead compensatcr is presented. 


Consider the uncompensated system shown in figure 
iJ. 
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4 
s(st+l) (s+2) 


Design Example 1, Block Diagram 


Figure IV-3 


The asymptctic Bede diagram for the open loop system is 
shown in ir figure IV-3A. As can te seen frem the cpen loop 
Bode plot, even thcugh the system is stable it is lightly 
camped. A single section lead compensator is to be used to 
increase tke phase margin. A complete discussion of the 
reasons for selection of the lead compensatcr and a detailed 
discussicn of this particular example is presented Ly Thaler 
and Brown [18]. fhe transfer functicn fcrm of the single 
section lead ccmpensator consists of a first crder numerator 
and a first order denominator. A desired magnitude and 
phase prcfile, consisting of ten discrete frequency foints 
over the range frem Q.2 to 20.0 radians, was selected to 
give a phase nargin cf approximately 45 degrees in crder to 
increase tke damping cf the system. As an initial guess the 
ccipensator was assumed to have a gain of 1.0, a zero at the 
Origin ard a ¢cole at ~1.0. This informaticn was suprlied to 
the prograr in the prescribed format and the optimized 
results in the form of magnitude and phase plots are shown 
in figure IV-3£. Here, as thrcughout the remainder of the 
example prcklems, the diasond shaped symbcls represent the 
desired values specified by the program user, the X symkols 
represent the values computed by the program at the discrete 
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frequency values specified using the compensator farameters 
returned as the sclution from the minimizaticn routine, and 
the continuots curve represents the response computed ky the 
prcgram again using the parameters returned as a solution 
but with a mcre dense set of frequency values used in crder 
to show any resonant peaks or anomalies in the magnitude and 
phase curves which might be difficult to visualize frcem a 
sparsely selected set of discrete values. As can be seen in 
figure IV-32 the results match the desired response well 
over the trarge cf interest. Figure IV-3C shows the same 
results using the magnitude versus phase plct to portray the 
open locg frequency response. The numerical values cf the 
ccmpensatcr fparameters returned from the program as a 
scluticn alcng with the starting values used are shcwn in 
meaqures J¥-3C, IV-sE, and IV-dF. A comparison with the 
scluticn values given in ref. 18 will show that they agree. 
Figure IV¥-3G shows plots of the differences between the 
actual ard desired open loop frequency response cf the 
resulting system at the discrete frequencies considered in 


tke minirizaticn process. 


This problem as descriked reguired aprroximately 500 
iterations sithirn the minimizaticn routine to achieve the 
sclution ¢hcewn. This however includes twe restarts to 
attempt insuring that a global minimum had keen achieved and 
that terminaticn was not due tc a local mininun. Frog the 
begining tc the first restart required approximately 250 
trials with the cost function being reduced from a value of 
6.23 to a value of 0.00711 =(recall that the minimun 
achievable urder ideal condition is zerc which due to 
machine eérrers alone in representing numbers would not 
likely ever ke achieved). The entire scoluticn requireé 32.5 
seccnds cf CEU tine, a good 10 seconds of which was used in 
the prodetcticr of Calcomp plots. These numbers ccncerning 
the time required for soluticn are presented for what they 
are wortk in that extrapolaticn of these numbers tc more 
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ccuplicated problems is not a straight fcrward task. As 
menticned in Chapt III the minimization routine does [frecome 
inefficient when the number of variable parameters is 
allewed tc kteccme large. It was felt that a detailed 
investigaticn of the efficiency and time required for 
sclution cf froblezs as a function of the number of variable 
parameters and the number of discrete frequency values tc be 
ceonsidered wae beycnd the score of the kasic developmert of 
the prograg. However, the amount of time required for 
sclution cf scme of the example problems that follow will be 
given in crder to give the user an "intuitive feel" fcr the 
approximate amcunt of time that may be need to _ kandle 


certain tyres cf problems. 


In crder to illustrate that the starting foint 
chcsen shculd have reasonable values, but that the exact 
values chosen should not effect the final results 
significartly, the same problem was run a second time with 
the pole and zero of the compensatcr bcth chcsen tc fe at 
the crigin. Again an assumed initial gain cf one was chcsen 
for the ccmpensator. The graphical results cf the solution 
are shown in figures IV-3H through IV-3J. The numerical 
results returned using these new starting values are not 
shewn here since they were identical to within three decimal 


places with those values shown previously. 
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2. Ccmrensator Lesign, Examp 


lo 
Ito 


{kis example problem (discussed in detail in ref.10) 
has been chcsen nct only to familiarize the reader with the 
use cf the program, but also to illustrate one of the 
shortccmings cf the design pregram. In particular, while 
the program was successful in finding the correct sclution, 
difficulties were encountered in the root finding 
subroutines. In this particular problem, during the final 
phases cf execution the roots cf a fourth crder pclynonmial 
must ke fcuré in order to calculate the phase response using 
suktroutine SIZFULA tvhich does these calculaticns using a more 
densely spaced set cf discrete frequencies than ir the 
pinimizaticn routine. Since numerical root finding routines 
have difficulty in finding roots £Or fourth crder 
Eclynomials, the lack of good convergence for the rocts of 
the system in Ecth rcot finding subroutines caused the phase 
calculaticns normally done by subroutine SIMULA tc be 
abcrted. It should be emphasized that failure of the root 
finding subrcutines in this case only effects the generation 
of the "ccntinuous" phase response and in ne way created any 
prcklems fer the sminimizaticn part of the prcegram as 
evidenced Ly the resulting correct values for the 
ccmpensatcr farameters. The original unccopensated systen 
being ccnsidered is shown in figure IV-4. A Bode diagram of 
the unccufpensated system (figure IV-4A) reveals that the 
system as it exists is unstable. If stability were the only 
consideraticr, a simple attenuation of the plant gain would 
be sufficient tc provide adeguate phase margin tc insure 
this. Hewever, as with many series compensation frotlems, 
here we must effectively strike a compromise between 
effectiveness cf ccntrol (generally asscciated with a higher 
gain) and stakility of the system (generally associated with 


a lower gain). As discussed in ref. 10a simple gain 
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attenuaticn swculd result in stability being achieved, but 
also has ar undesirable effect upen the velocity error 
constant. Tkus as discussed in detail in ref. 10, a lag 
compensation retwcrk appears as the best candidate to 


acccmuplish tke job. 






5°+20s°+100s 


Design Example 2, Block Diagram 


Figure IV-4 


In this particular problem twenty discrete frequency 
values over the range from 0.005 to 50.0 radians were 
selected ard the desired magnitude and phase values 
specified at these frequencies. As in the frevious eéxample 
proklem the form cf the compensator is selected as having a 
first crder numerator and first order denominator, that is 
in this case a sincle section lag compensatcr is to ke used. 
A type one fcrm of the cost function was selected for the 
Binimizaticn fprocess. The magnitude and phase profiles of 
the desired response along with the resulting response with 
the compensator included in the cpen loop system are shown 
in figures IV-4B, IV-4C, and IV-4D. As mentioned, the 
continuous curve for the fhase response was not cosputed 
because cf the failure of the reot finding recutines to find 
the cpen lccp system roots. It can be seen in figure IV-4C 


that the phase response does agree well at the discrete 
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frequencies specified. Since we are dealing with a mininuo 
phase system and figure IV-4B shows that there exists no 
resonant freaks in the magnitude curve the user may conclude 
that the phase response will also -te well tkehaved fcr _ the 
various frequency walues ketween the explicitly specified 
discrete values flotted. This case should SEIVE to 
illustrate that the program reguires that the user still 
prcvide ar interpretation of the results. The cogaputer 
output cf figures IV-4E, IV-4F, and IV-4G shows that the 
ccapensatcr farameters are in agreement with the freblen 
scluticn as resented in ref. 10 if the reader wishes to 
make the ccmparison. Figures Iv-4H and IV-4I show flots of 
the differences between the actual and desired magnitude and 
phase at tke discrete freguencies selected for the 
bBinimizaticn frocess. 
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Anctker lag type compensation proktlem is fresented 
in this examgle. This time, however, the existence of a 
fourth crder pclynomial for e¢ither the cpen loop systen 
transfer functior rumerator cr denominator has fpurfosely 
keen avcided tc insure that the program will run te 
ccipleticn. The uncompensated system blcck diagram is shown 
in figure IV-5. 










: FOO 
s (stl) 


Design Example 3, Block Diagram 


Figure IV-5 


A Bcde diagras for the uncompensated system along with a 
detailed discussicn of the problem may be fcund in ref. 17. 
The unccapfensated eystem does not possess an adequate fhase 
Bargin tc irsure froper system performance. Thus, a single 
secticn lag ccmpensator is to be used in ecrder to achieve 
the desired cpen loop frequency response. FCesired gain and 
Fhase prcefiles are selected for the ccmpensated ofen loorp 
frequency response using 20 discrete frequency points over 
the range frem 0.01 tou s00.0 Cadians. The prcfiles 
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selected alcng with the results achieved are shcwn in 
figures I¥-5A, IV-5B, and IV-5C. The numerical values 
assumed as initial estimates of the compensator parameters 
and the numerical results of the output are shown in figures 
IV-SD, IV-SE, and IV-5F. Plcts of the differences Letween 
the actual magnitude and phase, and the desired magnitude 
and phase are shown in figures IV-5G and IV-5H. In this 
particular prcklem approximately 1200 iterations of the 
girimizaticn routine and 1 minute and 25 seccnds were 
required for sclution of the _ problen. In order to 
illustrate the idea that for minimum fhase systems the 
ccmupensatcr parameters may be determined to within orly a 
fixed lcss cr gain when using only the phase profile, this 
same problem was solved again using the type 3 cost furction 
which ccnsiders only the difference between the desired and 
actual phase at tke discrete frequencies. The results are 
shcwn in figures IV-5I, IV-5J, and IV-5K. As can be seen 
frcm the magritude plot in figure IV-5I, there is a constant 
gain errcr in the resulting magnitude curve. The numerical 
values assumed at the start and those returned after the 
minimizaticn cf the cost function are shown in the couputer 
Output of figures IV~-5L, IV~-SM, and IV-5N. The resultant 
Magnitude response is higher than desired. The necessary 
cerrecticn tkat must be applied to obtain the desired 
response may ke read directly from figure IV-50, which shows 
the difference between the actual and desired magnitude 
response. That is, as can be seen from figure IV-5C, the 
Cresultant tagnitude of the ofen loop systen after 
ccupensaticn is approximately 36.8 db higher than what is 

desired. Since the program does not alter any of the plant 
Farameters this peans that the compensator gain value 
returned frcem the program when the type 3 ccst functioE was 
used must be decreased by approximately a factor of 69.2. A 
ccBfariscn of the compensator gain values shcewn in figures 
IV-SF and IV-5N will show that they due indeed differ by 
this value. Also the astute cbserver May nctice that there 


76 





is a slight difference in the compensatcr pole value 
returned as the sclution when the problem was solved using 
the type 3 cost furction. While the difference in valves is 
small it is suspected that in this particular case it is due 
to the fact that only a finite frequency range has been 
considered ard the lower frequencies have been "slighted" 
scuewhat in that the nonlinear phase curve has oct yet 
flattened cut at 0.1 radians. This is also the same 
situation tktat existed the first time the prcgram was 
executed, but in this case the magnitude prcfile was also 
taken into account in the cost function and any deviations 
cf the resultant magnitude values dominated the cost 
functicn value. In other words, it appears that ir this 
particular protlem the cost function is more sensitive to 
variations in the magnitude values than in the variaticns in 


the phase values. 
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Figure IV-5N 
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In this example a slightly more complicated problem is 
presented, in that the final soluticn requires a double 
section ccmpensator of the fcerm of a notch filter. A _ block 
diagram cf the unccmgensated system is shown in figure IV-6. 
As discussed Ey Thaler and Frown [18], from whom this 
example was taken, the uncompensated system is unstable and 
the Ecde diacram of the uncompensated system shown in figure 
IV-6A indicates a phase margin of approximately -30 decrees. 
With such a Ilarge negative phase margin £Cr the 
uncompensated system and the rapid drop in slcpe cf roth the 
gain and fghase curves the designer might suspect that 
ccmpensation using only a single secticn ccmfensator will be 
difficult to accomplish. Again desired magnitude and fhase 
prefiles are selected and initially a single section of 
ccmpensaticn is assumed to observe how close to the desired 
response tke system will perform. As can ke seen in the 
resultant craphical putput of figures IV-6A, IV-6B, and 
Iv-€cC, the single section compensator fails to meet the 
reguired frecuency specifications. In fact, while it is 
Ecssible to stabilize ths system with a single section 
compensatcr the bandwidth will be excessively large and the 
desired Bagnitude profile prevents the program fren 
acccmuplishing this. The numerical results returned frem the 
Single secticr ccmpensator run are shown in figures IV-6D, 
IV-6E, and IV-6F. As can be seen in figure IV-6F, the Routh 
test of tke characteristic equation also indicates systen 
instability which is to be expected. The magnitude and 
Phase difference curves of figures IV-6G and IV-6H, which 
indicate tte difference between the specified magnitude and 
phase prcfiles and the values actually achieved, indicate 
that additicnal compensation may be needed at the higher 
freguency tranges. Thus a double secticn compensatcr was 
assumed and the resudts achieved are illustrated in figures 
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IV-61, I¥-6J, and IV-6K. The program is able to satisfy the 
specificaticns with a double section of compensation and the 
parameters required to accomplish this are shown in figures 
IV-6L, IV-6¥%, and IV-6N. The differences betweer the 
desired and specified magnitude and phase values fcr the 
double secticn cf compensation inserted in the system are 
shcewn in figtres IV-60 and IV-6P. 
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this final compensator design example represents a 
prcblem cf a slightly different nature than those fresented 
to this feint. The servo system under consideration here 
bas a zero order held present in the error channel ard in 
addition, due to the large coefficient values that result in 
the plant transfer function, the problem requires frequency 
scaling in crder to insure that the numerical limitaticns of 
the particular computer that was being used would nct be 
exceeded. {he original single loop system is shown in the 
block diagrams cf figure IV-7. Series compensation of this 
system was desired in ocrder to meet the fcllowing 
specificaticns; 1.) the open loop magnitude at 377 
radians/seccrd 2 30db, 2.) the gain crosscver frequency 2 
18E5 radiansysecond, and 3.) the phase targin at gain 
crcssover 2 45 degrees. These specifications were to re met 
uSing only a series compensator and the parameters of that 
part of tte system shown in figure IV-7, other than the 
compensatcr, waS tc remain unchanged. An initial attempt at 
using the CALICO program to design a ccmpensator resulted in 
a warning of excessive number size during computaticr and 
reccmmended frequency scaling of the problem. Thus a scale 
factor of 10C0 was selected as a convenient value and the 
entire system, including the sampling frequency of the zero 
order hold, sas scaled down in frequency by the valve of 
this scale factor. Using Laplace transfcrm notaticn, this 
was acccaplished by a straight-forward substituticn of 
S = 1000S, where S represents the new scaled freguency 
values. A blecck diagram shcwing the values of the systen 
Parameters after scaling is given in figure IV-7A. When 
this type of scaling is performed the resultant numerical 
values cf the program need only be multiplied fry the 
appropriate ccenstant value in order to be corrected for use 
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in the unscaled systen. 


Criginally, work on this particular system indicated 
that a fifth order compensator could provide the necessary 
phase margin and gain crossover frequency, rut that the 30 
db specification at 377 radians/second was not satisfied. 
Thus, with this information as a starting point the 
Magnitude frofile cf the system with the fifth order 
ccmpensatcr in the loop was adjusted at the .low frequency 
end to satisfy the 30 db requirement at the scaled frequency 
value of C.377 radiansysecond. The phase prcfile was left 
unchanged. The results of the program can he seen in 
figures IV-7E through IV-7I. As can be seen frog the 
graphical cutput the desired specifications have teen 
satisfied Ly the resulting compensated systeg. 
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A. GENEBAL LISCUSSION 


In tke rocess of running various test problems to 
exercise and verify the algorithm presented in the preceding 
chapters, it was found that in additicn to its use as a 
ccmfensatcr design program, the basic method could also be 
applied tc the synthesis of transfer functions from 
freguency respcnse data. That is, if the magnitude and 
phase measurements of some system are cktained it is 
possible tc vse the program (CALICO) in order to oktain an 
approximate linear analytical expression for the transfer 
furcticn cf the system from which the measurements were 
obtained. This type of problem is quite cften encountered 
in. the mcdeling prcecess of dynamic systems. That is, the 
engineer may have measured input~output data of the system 
in the ferm cf a frequency response and is confronted with 
oktaining a transfer function which will mcdel this system. 
This may ke the result of the process being toc complex to 
be analyzed and the pertinent equations written on the basis 
of physical laws or the governing physical laws may nct be 
completely understood in the case of systems involving new 
technologies or reasearch efforts. Quite often the 
prccedure tcllowed in synthesizing or determining an 
approximate Jinear transfer functicn from frequency response 
data is te ccnstruct a Bode plot for the measured data and 
to fit the measured characteristic with straight line 


asymptotic approximations where possible and quadratic 
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factors where sharp peaks in the measured data occur. As in 
ccppensatcr design by classical frequency domain techniques, 
this precedure for synthesizing transfer functions can 
invclve a ccnsiderable amount of time consuming and tedious 
"cut and try" cn the part of the engineer. 


The algcrithm previously presented is so structured that 
the use cf tke program for this type of transfer furction 
synthesis is possible with no modifications to the coding 
itself. Since the algorithm was not originally developed 
for the specific purpose of transfer function synthesis 
there are several precautions that must be taken in crder to 
insure that the transfer function obtained does indeed 
adequately represent the system from which the measurements 
were oktained. In particular when using the program in 
crder tc synthesize a transfer functicn, waa sefarate 
Simulaticn frcgram for time domain response should be used 
to verify the transfer function obtained. MThis procedure of 
crecking the time domain response should be followed fcr any 
transfer function cbtained from frequency response data. 
The primary ccncern here is that it is possible to 
accurately represent the frequency response cf some systems, 
over a finite freguency range, by a transfer functicn of 
higher order than the true transfer functicn, but in some 
cases tke dcnrinant response to standard tiswe domain inputs 
will be ccnsiderably different than that of the actual 
systen. Alsc, while there is a direct relationshif between 
the time dcmein and freguency domain responses througk the 
inverse Fourier transforgo integral, this integration, 
strictly speaking, is taken over the entire frequency range 
frcm -~C tc +0, In practical applications it is ckviously 
necessary to limit the range of frequencies to be measured 
to some firite region. If, however, this regicn is not 
sufficiently large for a particular system the resulting 
transfer function obtained frem the frequency respcnse data 
Ray not acctrately describe the systen. fiwecnise 1S the 
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case, it should become readily apparent in a time dcmain 
simulaticn cf the equation obtained from the frequency 
response data. {The possible pitfalls in the use of this 
algcrithmw tc accomplish transfer function synthesis, such as 
the one just mentioned, will be discussed later when some 


example protlems will also be presented. 


A linear time invariant system may Le expressed as a 
ratio cf two frequency dependent pclyncmials of the forn 








D n= 
as +a Ss + n60 © das. fd N(s) 
n n-1 1 0 
H (s) me, nar — os 
bes) @ + 5 s + n56 °° DS tab CD(s) 
m n> 1 1 0 


where s is the frequency dependent ccmplex variable cf the 
Laplace transfcrm. Thus the value of this transfer function 
is dependert upon the value of s and the unknecwn 
coefficients cf the numerator and denominatcr. Note that 
this is analcgous to the situation that occurs in the Ccésign 
of a cascade ccmpensator for a unity feedback system with 
the plant represented by a unity gain block. Thus fcr the 
purposes cf using the program, if we imagine the measured 
frequency response data to represent some desired system's 
open loop frequency response profile, and represent the 
Flant of that system by a constant gain of unity, then the 
transfer function H(s) to ke synthesized may ke considered 
the same as a cascade compensatcr in the system. Thus by 
representing the problem in this manner the CALICO program 
may be used to approximate the transfer function cf the 
actual systesz Fy varying the parameters of this pseudo 
ccmpensator (actually the transfer functicn of the systen) 
in crder tc sinimize the cost function. The cost function, 
bere represents the error Letween the measured frequency 
response data (entered as the desired refile) and the 


freguency resfense simulation calculated by the pfrcgran. 
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Thus in effect the program minimizes the errcr between the 
transfer function approximating the system and the measured 
freguency respcnse of the system at the discrete freguency 
points which form the profile infut into the progran. 


As is the case when using the program tc design 
ccupensatcrs the user must assume a specific order fcr the 
numerator and dencminator representing, in this case, the 
transfer function to be synthesized. While unfortunately 
the program will not vary the order of the numerator and 
denopinatcr pclynomials in order to obtain the best match to 
the measured frequency response, the user can employ the 
measured Gata in order to make an initial guess, as it were, 
of the form cf the transfer function. After this, one may 
scrutinize the results from the program, which grafts the 
measured and simulated responses, to make a decision if the 
fcrm of tke transfer function should be adjusted in crder to 
more closely represent the system. Here alsc, a time dcmain 
Sigulaticn sight ke extrerely useful in guiding decisicns 
regarding the particular form cr modification that should he 
tried in crder to obtain a "good" model for the system under 
investigaticn. The fact that the user selects the specific 
crders cf the transfer function numerator and dencminator 
and the minisizaticn routine varies the parameters to cbtain 
a BSinimup ccest function while leaving the crders unchanged 
presents the potential of obtaining a low crder model of a 
higher crder system. While for linear systems, with known 
transfer functions, there are numerous analytical techniques 
available for accoaplishing this, the capability exists here 
of essentially oktaining a reduced order transfer function 
without a friori knowledge of the higher crder transfer 
function describing the system. The only data required by 
the program is the frequency response of the higher order 
System and ncthing concerning the parameters of the high 
crder transfer function is required. 
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Since scme transfer functions may contain roots of the 
nuserater feclyncmial in the right half of the s~plane 
(ncnminigsum phase systems) this type of transfer function 
should ect ke overlooked when attempting tc model a systen. 
fo this erd a nonminimum phase cption may be chosen ir the 
CALICO cceputer algorithm. When this opticn is selected the 
search area fer the numerator coefficients is expanded to 
include Eoth positjve and negative values and the implicit 
censtraints requiring left half plane rcots LOE the 
nuperator feclyncmial are eliminated. Selection of this 
cption dces nct guarantee that a nonminimum fphase transfer 
function will result from the CALICO program but it does 
allow this type of transfer function to be considered as a 
feasible sclution cf the problen. 


B. TRANSFER FUNCIION SYNTHESIS: EXAMELE PRCELEMS 


In this section several example prceblems are presented 
to illustrate the use of the CALICO program in determining 
transfer functions from frequency response data. These 
€xamples illustrate the validity of the kasic concert of 
using the minisization algorithm to determine the transfer 
furction, tut they are by no means an exharstive 
Presentation of all possible variations of cost functions 
and freguency ranges that say ke employed in determining a 
transfer function representation of a system. For example, 
as will ke shcwn in the example problems, different transfer 
furction representations may be obtained ky considering 
different ranges of the measured frequency response. That 
is, while the order pf the numerator and denominator are 
fixed by tke user, the parameter values returned frcem the 
prcgram ray vary depending upon the frequency range being 
considered. This may be desirable if seperate low frecuency 
and high frecuency representations are being determined. 


iz 





Also in using the type one cost function, which includes 
both phase and magnitude, the cost surface is different than 
if the tagritude difference alone were used in the cost 
function. This has an effect on the resulting farameter 
values, particularly when lower order transfer functions are 


chcsen tc represent higher order systems. 


The necessary data and the input format employed are 
almcst identical to those used when the prcegram functicns in 
its originally designed compensator design mode. Figure V-1 
illustrates the input data deck used for the first synthesis 
exanuple problem. As can be seen the only significant change 
is that the slant has been represented by a unity gain 
transfer function. This has the effect that the resulting 
"ccmpensatcr" represents the transfer functicn modelinc the 


entire syste. 


> 
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1. Synthesis Example 1. 


2 ee oe 


fhe freguency response data, shcewn graphically in 
figure V-Z, was generated from a known transfer function 
using a ccmputer frequency response program. As can he 
seen, twerty magnitude and phase values were selected at 
discrete freguencies over the range from 0.1 radians/second 
to 30 radiansysecond. Examination of the frequency response 
data suggests a transfer function of the form of a corstant 
gain over a g¢clyncmial of some unknown order. This foro 
will yield a constant magnitude at low frequencies anda 
decreasing wagnitude as the frequency increases. Clcoser 
exaBinaticn cf the magnitude curve in the range between 2 
and 20 radiansysecond indicates a slope of approximately 60 
Gkyéecade, igplying a polynomial of at least third crder in 
tke denoninatcr. The phase curve extends frem 10 degrees to 
330 degrees over the range cf frequencies ccnsidered. This 
indicates that the denominator may be of fcurth crder or 
higher. As a first guess a third order denominatcr is 
assumed. The results from the CALICO program are shown in 
figure V-2A through V-2H. The initial guess for the 
denominatcr rects and the system gain were chosen as 1, as 
Can be seen in figure ¥-2F. Figures V-2A and V-2D show that 
the resulting magnitude curve approximates the measured 
response well at the lower frequencies with a slight 
deviaticn becining tc appear at the higher frequencies. The 
phase plct cf figure V-2B and the plot cf the difference 
between the actual phase and desired phase shcwn in figure 
V-ZE indicate a significant amcunt of deviation frcm the 
measured values in the higher freguency region. This 
Suagests that ancther pole in the transfer functicn is 
necessary tc mcre accurately model the systen. Thus, the 
prcegram ‘as again executed, this time using a fourth order 
dercminatcr. The results illustrated in figures V~2I 


i322 


thrcugh V-2EF indicate that this second form of the transfer 
functicn gcre accurately represents the systen. Since the 
Criginal transfer function is known in this case, an 
opcrtunity exists to evaluate the performance of the prograr 
in terms cf the accuracy of the transfer function returned 
after minimization of the cost function based on the input 
frequency respense data. This of course would generally not 
be the case when trying to determine the transfer furction 
of an actual system that was to te modeled. fhe parameter 
values returned frcem the CALICO program shown in figure V¥-20 
may be ccmpared with the original transfer function given 
below. 
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t(s) = ——_-—_- -—_-_—_ --—__—_ + -— 
(s + 1) (s + 2)(s + 5) (s + 10) 


As can Ee seen, the values of the transfer functicn roots 
returned frcem the pregram for the fourth order denominator 
are in agreement with the values of the original transfer 
function used to generate the frequency resgcense data. Lf 
for the system under consideration only the magnitude 
response were cf interest, the transfer function ccnsisting 
of a ccnstant gain over a third crder dencminator nicht be 
of sufficient accuracy to model the system. This judgement, 
of course, rests with the user who is familiar with the 
intended use of the transfer function modeling the system. 
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Figure V-2P 





As a second synthesis example, consider the 
freguency response data shown in figure V-3. Here the 
magnitude and phase response for a system are shotn at 
twenty distinct frequency values over the range frcem 0.1 tc 
100C.0 radiansyseccnd. As can be seen from this figure the 
positive slcpe of 20 db/decade of the mgnitude curve at the 
lower freguencies suggests a transfer functicn with at least 
one zere. This is also confirmed by the fhase resfense, 
which at the lower frequencies initially has a value cf 90 
degrees. At the higher frequencies inspection of the 
magnitude plot shows a negative slope of apffroximately 40 
dktydecade, indicating the dominance of a third order 
dencminatcr in this area. The phase plct approaching -180 
degrees at these higher frequencies alsc suggests this. 
Thus, initially a form for the transfer function consisting 
of a first crder numerator and a third order denominatcr was 
selected tc sodel the system from which the frequency 
tesponse data was obtained. The results are shcwn in 
figures V-3A thrcugh V-3H. As can be seen fron the 
graphical ottfut presented in the figures, this particular 
form of transfer function appears to match the wneasured 
freguency response guite well over the entire range of 
frequencies ccnsidered. The rcots and gain cf the transfer 
function rettrned from the program are shown in figures V-3G 
and V-3H. As in the previous example it hapgrens that the 
€xact transfer function used to generate the measured data 
input tc tke cCALICC program is known. This transfer 
function is given Felow in order that the reader may ccmfare 
the values cf the orjginal transfer function parameters with 
these returned frcm the progran. 
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9000s 


Tt(s) = ———__________ 
(s + 2.4) (s + 5) (s + 150) 


As can be seen fre the data presented in figures V-3G and 
V-2H, the parameter values returned from the program are in 
clese agreement with the original transfer function values. 
While it is urlikely that they will ever be exactly the same 
as the criginal values there is certainly nct a significant 
difference between the two sets of values fcr most 
engineering furroses. 
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Figure V—3C 
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3. Synthesis Example 3, 


A frequency response of a system with resonant 
peaks, shcewr in figure V-4, is ccnsidered in this example. 
The same example is considered in ref. 19, where a 
ceombinaticn cf graphical and analytical techniques are used 
to deterpine the transfer functicn. Initially a transfer 
function ccnsisting of a first order numerator and second 
order dencuinator was assumed in trying to mcdel the system. 
This chcice was based on the low and high frequency 
Magnitude curve slopes of approximately +20 db/decade and 
-20 db/decade respectively and the low and high frequency 
phase values cf approximately +90 degrees and -90 decrees. 
This fors failed to match the dip in the magnitude curve in 
the regicn cf 10 radians/second and the jump in the phase 
curve in the same frequency area. After increasing the 
crder of the numerator and denominator, a transfer furction 
consisting cf a third order numeratcr and fourth order 
dencminatcr was found to provide a phase and gain curve that 
closely approximates the measured frequency response. The 
results cf this are shown in figures V-4A through V-4H. The 
nuserical values of the transfer function pararéeters 
returned frc@ the program clesely approximate those of the 
actual system from which the frequency resfonse measurements 
were obtained. The interested reader may verify this fact 
by comparing the values shown in figures V-4G and V-4E with 
these of the actual system given in ref. 19. 
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In this final example problea a low crder 
approximaticn to a known higher order system is illustrated. 
The original system transfer function under consideration 
consists of a seventh order numerator and an eight order 
dencminatcr given by 


7 6 s & 3 2 
16s +4E3s +6010s +36380s +122664s 4222088s +185760s+40320 


8 ? 6 s 4 3 2 
Ss +36s +546s +4536s +22449s 4+67284s #+118124s +1095845+40320 


An analytical method of approximating this high crder 
transfer function is discussed in ref. 21. Initially 40 
discrete fpcecints, representing the frequency response cf the 
high order transfer function were selected OVEL the 
freguency Lange £rom 0.01 radians/second tc 100.0 
Tadians/second. In order to be able tc form a ccmparison 
with the reduced order model obtained in ref. 21, a transfer 
functicn ccnsisting of a first order numerator and a_ second 
crder dercninator was assumed. The results returned from 
the program tsing this form of transfer function are shown 
in figures V-5A through V-SH. It is not surprising that 
there is a significant amount of deviation in the sagnitude 
curve in the vicinity of the resonant effects of the higher 
order system nor that the phase response cf this’ reduced 
order mcdel dces not match the higher order phase response 
in the hicher frequency range. In order tc compare these 
results with analytical methods for forming reduced order 
ucdels, the frequency response for the reduced order snodel 
of ref. 2z1 is shown in figures V-5I through V-5M. The 
nugerical values returned from the CALICO prcgram are cf the 
same crder cf magnitude as those found Ly the analytical 


technigue, rut as can be seen frem the magnitude curves of 
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the two apprcaches (figuresV-5A and V-51I), the analytical 
method matches the higher order transfer function extremely 
well at the lower and higher frequencies and makes little 
attempt tc match the higher crder system in the rance of 
freguencies at which the rescnant peaks occur. The results 
returned frcemp the program, on the other hand because they 
are a functicn cf the deviaticn of the response frcm the 
desired values over the entire range of frequencies, tend to 
distribute any error between the responses over the entire 
range of frequencies being considered. A slightly different 
representaticn of the system may ke oktained by varyirg the 
frequency rarge under consideration and the type of cost 
function used in the program. For example, figures V-5N 
through V¥-5U illustrate the effects of using a type twe cost 
furncticn and limiting the frequency range from 0.1 to 100.0 
tadians/seccnd. Again the decision as to which acdel is 
best must te made by the user in light of the intended use 


of the recuced order representation. 
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Figure V-5H 
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Computer Numerical Output, Run #2 


Figure V-5U 





The investigations carried out in this thesis illustrate 
the feasibility cf automating the classical design of 
compensatcrs by specifying a set of desired values fcr the 
magnitude and fhase response of the open lcop system at a 
nugker cf discrete frequency values over the range of 
interest. The results achieved show that the ccmplex 
girimizaticn methcd of M. J. Box is quite capatle of 
handling the ccnstrained Minimization problen and 
conveniently provides a means of avoiding the numerical 
differentiation cf the cost functicn, required with gradient 
search technigues. The incorporation cf staktility criterion 
in the inplicit constraints of the minimization algcrithn 
alsc insures that the program will not design an unstable 
compensator. Ey specifying the desired respfense in terms of 
the magnitude and phase profiles, ccmmon frequency domain 
specificaticns such as open loop bandwidth, phase margin, 
and gain margin have been incorporated into one ccmmon 
format and tke necessity for specialized routines to check 
these specificaticns have been avoided. Those example 
prcklems presented for which known results had been 
estaklished Ey other techniques show that the algcrithn 
presented prcvides accuracy that is certainly within the 
necessary tclerances for most engineering applicaticns of 
compensatcr design. 


Also, while not specifically designed for the purpose, 
the basic fremise cf the algorithm's Binimizing the 
difference tetween some desired frequency response and the 
frequency respense of the system as calculated by the 
ccaputer prcgram while the transfer functicn parameters are 
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varied has been shcwn to be effective in both determining a 
transfer functicn from measured freguency response data and 
in forming lcw order models of known higher order linear 
systems. In the formulation of lcw order mcedels the choice 
of cost functicn and the frequency range specified may he 
usec tc emphasize different aspects of the frequency 
response which tke low order model is tc approximate. 


While this werk has shown that the kasic algcrithn 
presented can be used to effectively design series 
cc@fensatcrs tased on frequency response specifications, 
consideratle work remains to be dcne in this area. 
Specifically, an option to select that only real foles and 
zercs be ccrsidered as viable solutions to the compensator 
design algorithm would not be a difficult additicn tc the 
program and should prove quite useful. Also, the 
inccrporaticn of a provision for handling pure time delays 
would not reguire significant modification cf the progran 
ceding and wculd extend the capability of the algcrithm tc 
the soluticn of problems involving transport delays. 
Investigaticn of the method presented should be further 
pursued in terms of multiloof systems and systems invclving 
ccmpensaticn in the feedback loop. It may also prove 
worthwhile tc investigate extending this methcd to nonlinear 
ccmponents using describing functicn technigues. 


In the areas of transfer function synthesis and reduced 
order mcdeling, a more detailed investigaticn of both the 
effects cf different cost functions and considering 
different frequency ranges Should be conducted. In 
addition, it appears feasible to weight the cost furction 
over a particular frequency range by varying the density of 
the discrete freguencies at which the desired magnitude and 
phase prcfiles are specified. Time did not fpergit a 
detailed investigation of this idea and further 


investigatior intc this area may be aprlicable to the 
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program utilization in both the areas of compensator 
and transfer function synthesis. 
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APPENDIX A 


DESCRIPTION OF SUBROUTINES 


In this section a brief description and backgrourd of 
various subroutines is presented to aid the interested 
reader in better understanding the internal logic of the 
prcgram. Also these descriptions are intended to facilitate 
the. task of anycne desiring to make mecdificatiors ofr 
additions tc the program in order to handle systems cf a 
configuraticn cther than the unity feedback form assumed 
thrcughout this work. 
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SUBROUTINE CONORM 


The fpurgese cf this subroutine sub-pfrogram is to 
Nnorwalize tte coefficient cf the highest order term to 
unity. Mbis is accomplished in a simple, straight fcrward 
manner by dividing all the coefficients of the polynomial by 
the coefficient cf the highest order term. Thus a 
EClynomial given by 





n 
1 
P(s) = a es 
a: 
i1=0 
beccmes cn return from this subrcutine 
n 
a, e 
1 2 
P(s) = ) es 
a 
. n 
1=0 


n 
where a is the coefficient cf the s term. The normalizing 
n 


factor a is also returned to the calling prcgran. 
D 


It shcvuld be noted that the coefficients of the 
norgalized felynonial are returned under the same variable 
nawe as tke input coefficients. 


Definiticns cf the input and output parameters are: 


Input Variables 


es -————— A real one-dimensional array containing the 
ccefficients of the polynomial whose highest crder 
ccefficient is to be normalized to unity. JUpon 
return from the subroutine this array will ccntain 


the ncrmalized coefficients of the felynomial. 
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Ho o--~-------- The dimension cf the array containinc the 


ccefficients tc be normalized. 


VNOEM ------- The normalizing factor, returned tc the 
calling prcgram, by which all coefficients have keen 
divided. 
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SUBROUTINE PHASE 


This sukroutine sub-program is used to calculate the 
Fhase angle cf a pclynomial at a particular frequency. In 
order tec circumvent the difficulty associated with standard 
arctangent rcutines, which return values limited tc _ the 
range —"i£90 2, the phase contributions from the 
individual factors are calculated and then summed to form 
the total phase angle at a particular freguency. That is, 
if the pclyncnial is factored and represented as the product 
eof first crder terzs of the form 


n 


P(s) = Tl (s +z.) 
1 


i=0 


where in general z may be a complex number of the form 
i 
Z=a +gjk, then at a particular frequency s = je the 
i a i 


phase contrikution of the individual factors may be 


evaluated as 


(bo tw ) 
x 


a e 
as 
and the tctal phase angle at the particular frequency is the 
¢um of these individual phase contributions. The phase 
contributicn cf each of these first crder terms will be 
within the principal values of the arctangent routine cf the 
ccaputer and the problems generally associated with 
determinaticn cf the phase for higher than first order terms 
are avoided. 


Definiticns of the input and output parameters are; 


Input Variables 


, were er---- A reak one-dimensional array containins 


real farts cf the rocts of the polyncspial. 


a A real one-dimensional array containing 


isaginary parts of the roots of the felynomial. 


mene eo -- An integer value specifying the number 


rects cf the polynomial. 
oA -—---—~— The frequency in radians/sec at whick 
value cf the phase is to be computed. 
Output Variables 


oo --——-— The total phase contrikution of 


pclyncnial civen in degrees. 
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SUBROUTINE PINVRAT 


The fpurgpese cf this sukroutine suk-fregram is to 
rearrange tke coefficients of a polynomial in inverse order 
from that fed into the subroutine. Consider for example a 
Eeolynomial arranged in descending powers of s 

E(s) = a =. + a oo +t sss t+ aS £ 2 
D n-1 i 0 
Upen return frcm this subroutine the coefficients wculd be 
arranged in ascending powers cf s. That is the pclyromial 
weuld be arranged as fcllows 
E(s) =a +¢+¢a4S + «22 + 4 gr + a Ss 
G i n-1 n 
While kcth representations are mathematically equivalent 
the sequence cf the coefficients must be taken into acccunt 
when using various pclynomial manipulation techniques tsithin 
the computer. The method used in reversing the sequence of 
the input ccefficients consists of temporarily storing the 
first coefficient and sequentially shifting the remaining 
coefficients. The coefficient in temporary storage is then 
Flaced in the last coefficient location. This process is 
continued with the use of nested do loops and each time the 
shifting operation occurs the coefficient in temporary 
storage is laced one location lower than the previcus 
number retrieved from temporary storage. 
\ 

This stkroutine requires no cther subroutines or 
functions tc carray out the procedure. The size of the one 
dizensicnal array of coefficients that may te rearranged is 
determined ky the dimension cf the array in the calling 
progran. 


Definiticns of the input and output parameters are: 
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Input Variables 


a ----- ~~ A real one dimensional array containing the 
ccefficients tc be rearranged. Upon return frer the 
subroutine this array contains the input coefficients 


in reverse seguence from that entered. 


N er erer rr The dimensicn of the array containing the 


ccefficients tc be rearranged. 
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SUBROUTINE PMPY 


The purpose of this subroutine suk-progragz is to 
multiply the ccefficients of one polynomial X(s) times the 
coefficierts cf anether polynomial Y(s) to obtain the 
coefficients cf the resulting polynomial Z(s). That is, if 
it is desired to multiply two polynomials of the form 

2 n-1 n 


an * a Set a SS” tft aee t a s + as 
0 2 2 Toa D 


P (s) 


and 


© al n 
bo +¢+bst+bs * «ee + Dd s + kos 
0 a € n-1 n 


P(s) 


tcgether in crder to obtain a third resultant polynomial of 
the general fcrm 
P (§) =c *#cC¢S#C a + 56 t C a + Cc a 

3 0 1 2 i-1 i 

this sctkrottine is used in acccmplishing this. The order 
of the resulting pclynomial is calculated by forming the sum 
of the oréers cf the two input pclynomials as 

l1=no +n 

The varicus ccefficients of the resulting pcelynomial, which 
are returned ir ascending fpowers of the independent 
variable, are calculated as the sums of products of the 
input ccefficients by means of two nested do-loops. The 
arrangement cf the two do-loors is such that the approfriate 
product terms are summed in forming the respective 
coefficierts cf the resultant pclynomial. All input and 
output ccefficients must be real numbers. The size cf the 
Eclynomials that say be multiplied are contrclled by their 
respective dimension statements within the calling prcgraa. 


Definiticns of the input and output parameters are: 
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Input Variables 


) i re A real one dimensional array containing the 


ccefficients of the first polynomial tc be 
multiplied. These must be arranged in ascending 
powers of the independent variable. 


------- The dimension cf the array X c ontaining the 


ccefficients of the first polynomial. 


oo s---- A real one dimensional array containing the 


ccefficients of the second polynomial tc be 
multiplied. These must be arranged in ascending 
pewers of the independent variable. 


eS oss The dimension of the array Y containing the 


ccefficients of the seccnd polynomial. 


Output Variables 


------- A real one dimensional array containirg the 


ccefficients cf the resultant pclynomial. These 
ccefficients are arranged in ascending powers cf the 


independent variable. 


------- The calculated dimension cf the array Z 


ccntaining the coefficients Of the resultant 
pclyncpial. 
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SUBROUTINE PVAL 


The purfpese cf this sukroutine sub-frcgram is to 
evaluate a pelyncnzial having the complex variable s = o + jw 


and real coefficients a. For a given polyncmial of order n 
2 


of the fcrao 


n 


1 
P(s) = ) aes 
al 


i=0 


s is declared as a complex variable within the subrcutine 
and the indicated operations of summation and multiplication 
are carried cut by a simple lcop operation taking advantage 
cf the asscciative and distributive laws cf multiplication 
and additicn. This subroutine requires nc other 


suk-programs tc perform its function. 


Definiticns of the input and output parameters are: 


Input Variables 


oo ———-—--—-—— A veal one-dimensional array containing the 
ccefficients cf the polynomial that is te be 
evaluated. The coefficients are assumed to be 
arranged in ascending order according to the fowers 
of s with zero explicitly input for the coefficient 
of any missing powers of s occuring in the 
pclyncrial. The maximum size of this array is 
ccntrclled by the dimension of the array withir the 
Main frogram. The values of the coefficients in this 
array remain unchanged after execution cf the 
sukroutine. 
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a ————--=-- The order of the polynomial tc be evaluated. 
This value is one less than the dimensicn of A(I). 


a —————— The value of the real fart cf the ccnuplex 
variatle s at which the polyncnwial is to be 
evaluated. 

a —— — The value of the imaginary part cf the 
ccuplex variable s at which the polynomial is to be 
evaluated. 


Output Variables 


== ~~ The real part of the complex resultant value 
of the polynomial evaluated at s = o + jw, 


—————————-~——= The imaginary part of the ccmplex resultant 
value cf the polynomial evaluated at s =o + jo. 
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SUBROOTINE ROUTH 


fhe furgpcese cf this subroutine sub-frogram is to 
determine if a given input pclynomial possesses right half 
Flane rocts ty using the Routh stability critericn. The 
EClynomial ccefficients are input in descending pewers of 
the independent variable. From these coefficients the Routh 
array is calculated and the first column of this array is 
searched for an indication of the existence of fright haif 
Flane rocts. Thus from an input polynomial cf the forn 


n n-1 
E(s) = bs +b Ss + cee + bStd 
n n-1 1 0 
tke Reuth array 
st 
1 row b b L Stee os 
n n-2 n-4 
nd 
2 Low b b b siecaneiere 
n-1 n-3 n-5 
rd 
3 Low a a a sierecsuaiere 
31 32 33 
th 
4 row a a a eueveieiieiee 
41 42 43 
th . ; i 
n row a a a aateuel ciate 
n1 n2 n3 


st 
(n+ 1) rcw a a a 
(n+1)1 (n+1) 2 (n+1) 3 
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is formed, where in general the entries from the third row 


cnwerd are given by 
ie (n> 1} il? (n-2) (i+1)] / [2 (n-2) il? (n-1) (i+) 


a 
(n-1) 1 
The first cclumn of this array is then tested to determine 





Es 


if any of tke elesents are less than zero. Note that since 
all coefficients are limited to positive values by the 
search area restrictions in cther parts of the prcegram no 
prcvision has been made to handle the case cf the 
coefficients cf the polynomial being all negative. Thus a 
EcClynomial sith all negative coefficients should be 
aultiplied ky -1.0 prior to using this particular 
subroutine fer testing for right half flane roots. 


Definiticns of input and output parameters are: 


Input Variables 


Y 2---------- A real one dimensional array containing the 
ccefficients of the polynomial in descending fewers 
of the independent variable. 


N ----------- The dimension of the array Y. 


Output Variables 


Ay >----— Integer output indicating the absence or 
presence of right half plane roots 
O-~-Indicates the polynomial has no right half 
Flane rpots. 
1---Indicates the polynomial has fright half 
Flane roots. 
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SUBRCUTINE SEMEL 


The furpese cf this subroutine sub-frogram is to 
calculate the ccefficients of a polynomial from the reccts of 
the polyncaial. The calculation is done using ccmplex 
arithmetic in order to te able to compute the polyncmial 
coefficients when complex roots are present. The resulting 
FClyncmial ccefficients are, however assumed to he real. 
Thus, if ccmplex roots are present they must be in conjugate 
fairs. Alec the coefficient of the highest order term is 
set egual to unity. Consider a polynomial of degree n 
having n rocts. This gay be expressed as 


n 
P(s) = | | (s + 2.) 
1 

1=0 


where in general the z's may be complex. The coefficients 
1 


of the polyncmial that results when the above indicated 
multiplicaticn is performed may be computed Ly the fcruaula 


n71 
A. = (~1) e ) (product of roots taken n-i at a time) 
i 


Thus the coefficients of the pclynomial are given by 
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a = Z2 @Z @Z @®Z #®...°2Z 
1 2 3 4 


n 
a = 2 @2 #2 #22292 + Z @Z @®Z @©...°2Z 
1 1 2 3 nest 2 3 4 n 
a = Zz £2 + Z + -Z vez, Te 616) fe 
a = 1 
D 


The largest crdeyr polynomial that this particular sukrcutine 
th 
will handle is a 20 order. If larger pclynomials are to 
be considereé the dimension cf R within the subroutine must 
be changeca. 
Definiticns of the input and output parameters are: 
Input Variables 
a - -———— A real one-dimensional array containing the 
real farts of the polynecnmial rocts. 


a ————— A real one-dimensional array containing the 
imaginary farts of the polynomial rcccts. 


ee An integer specifying the order of the 
pclyncpial or equivalently the dimension of the 
arrays containing the real and imaginary parts cf the 
rcects. 


Output Variables 


——— A real one-dimensional array containing the 


resulting ccefficients cf the polyncgrial arranged in 
ascending powers of the independent variatle. The 
dinension of this array must be one larger than the 
digensicns of the arrays containing the real and 
imaginary parts of the roots. 
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APPENDIX B 


PROGRAM LISTING 
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APPENDIX C 


CESCEIPTION OF FORTRAN VARIABLE NAMES 
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A (I) 


ADAMG (I) 


AMAG (I) 


ANGDEL 


ANGDEE 


ANGNUL 


ANGNUE 


E (IL) 


CO (I) 


COF (I) 


COFASN (I) 


UNIT 


None 


db 


db 


LTegrees 


Radians 


Tegrees 


Radians 


None 


None 


None 


None 


DESCRIPTION 


Arra ee the coefficients 
of he plant ransfer function 
numerator . 


Array containing the desired cfen 
loop magnitude respense in db at 
the discrete frequencies specified 
by the designer. 


Output array containing the 
magnitude of the cren Loop 
frequency response in db at tte 
discrete frequencies specified. 


Phase angle contribution of the. 
open loop system transfer function 
denominator at a specific 
frequency. 


Phase angle contribution of tke. 
open Tooe system transfer function 
denominator at a specific 
frequency. 


Phase angle contribution of the. 
open loop system transfer function 
numerator at a specific frequency. 


Phase angle ccntribution of the. 
open loop system transfer function 
numerator at a specific frequency. 


aca containing the coefficients 
of the plant transfer function 
denominator. 


One dimensional working array used 
primarily in reading ee 
coefficients and infuting data to 
root finding subroutines. 


One dimensional werking array used 
Primarily to store resulting 
coefficients used in checking 
accuracy of root finding 
subroutines. 


ae Lay containing the coefficients 
of the cpen loop transfer function 
numerator. 
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FORT 
VARIA 


toto 
4 ts 
ral oe 


COFASL (I) 


COFWRK (I) 


DELPHI 


DELPSI 


DIFF 1 (I) 


DIFF2 (1) 


DMAG (T) 


DCPHSC (I) 


DWFQ 


FREQ (I) 


GAIN 


GAINCC 


UNIT 


None 


None 


Degrees 


Degrees 


db 


Tegrees 


None 


Degrees 


Radians 


None 


None 


None 


DESCRIPTION 


Array containing the coefficients 
of the cpen loop system transfer 
function denominatcr. 


One dimensional werking array used 
primarily to temporarily stcre 
poesee te coefficients in 
performing polynomial 
manipulations. 


Difference between two numerator 
phase angle ccntrikutions at. 
Successive discrete frequencies. 


Difference between two denominator 
phase angle contrikutions at. 
successive discrete frequencies. 


Difference between the desired and 
actual magnitude of the open OG 
frequency response at the discrete 
fhe quency values specified in the 
input. 


Difference between the desired and 
actual phase of the open lecp 
frequency response at the discrete 
fred tency values specified in the 
input. 


Array containing the desired 
magnitude response at the discrete 
frequencies specified by the 
designer. 


Input array containing the desired 
phase profile cf the open loor . 
System at the discrete frequefrcies 
chosen. 


Interval to be used if discrete 
Se CaueNCY values are to be 
incremented linearly. 


Logarithm of the discrete ; 
frequencies used in constructing 
abscissa for Bode frlots. 


Gain of plant transfer functicn. 


Gain of compensatcr transfer 
function. 
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IDB 


IERBXE 


INFORSB 


IPLOT 


IZOH 


KNOW 


NEODE 


NCD 


NCN 


NICHCL 


NMINFS 


NOMEG 


NSD 


UNIT 


None 


None 


None 


None 


None 


None 


None 


None 


None 


None 


None 


None 


None 


DESCRIPTION 


Input integer to be set egual to 
unity if the desired open lcor 
Magnitude response is input in db. 


Output integer indicating the 
completion code of the mininization 
routine. 


Alphabetic input indicating if 
coefficients are te ke read ir 
factored cr polyncmial forn. 


Integer input, to ke set to unity 
2£ CALCOMP plots of results are 
desired. 


Integer input set egual to aD gd 
there 1s a zero ordér hold in the 
error Slane path prior to the 
plant and compensatcr. 


Integer input indicating if 
discrete frequency values are to be 
read in or incremented linearly 
from the minimum freguency. 


Integer input to Fe set equal to 
unity if Bode plot is not desired. 


Order of compensatcr transfer 
function denominatcr. 


Order of compensatcr transfer 
function numerator. 


Integer input to ke set equal tc 
unity if Nichol's flot is not 
desired. 


Integer input to ke set equal to 
ee ape a nopaininum phase 
solution is to be allowed. 


Integer cc ead the total number 

of discrete frequency points teing 

considered in the infu eso eased 

be a 800 profile. This numEer must 
e . 


Order fo the system open lccrp 
transfer function denominatcr. 
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NTFD 


NTIFN 


NYQST 


N3 


N4 


CDENT 


CDENR 


CMAG (TI) 


CNUMI 


CNUMR 


CPHSD (1) 


UNIT 


None 


None 


None 


None 


None 


None 


None 


None 


Ncne 


None 


None 


LTegrees 


DESCRIPTION 


Order of the system open lccp 
transfer function numerator. 


Order of plant transfer functicn 
denominator. 


Order of plant transfer functicn 
numerator. 


Integer input to ke set equal tc 
unity if Nyquist pict is not 
desired. 


Dimensicn of array containing plant 
transfer function numerator 
coefficients. 


Dimensicn of array containginc 
plant transfer function denominator 
coefficients. 


Value of the immaginary part cf the 
systen epee aoe Lanster functicn 
denominator evaluated at a 
particular frequency. 


Value of the real fart of the |. 
systen open Scan transfer function 
denominator evaluated at a 
particular frequency. 


Output array containing the 
Magnitude of the frequency resfense 
computed for the oren Peete Ste en 
transfer function at the discrete 
frequencies specified. 


Value of the Meas naLy part cf the 
System open loop transfer function 
numerator evaluated at a particular 
frequency. 


Value of the real part of the . 
System open loop transfer function 
Numerator evaluated at a particular 
frequency. 


Output array containing the phase 
of the frequency response of the 
open loop System at the discrete 
frequencies specified. 
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CPHSR (1) 


SAMAG (I) 


SAPHST (I) 


WFREq(¢1I) 


WMAX 


WHIN 


WS 


XS (TI) 


UNIT 


Radians 


db 


Degrees 


Secs 


Radians 


Radians 


Radians 


Radians 


None 


DESCRIFTION 


Output array containing the phase 
of the freguency response of the 
open loop system transfer function 
at the discrete frequencies 
specified. 


Array containing the magnitude of 
the Simulation of the ES. 
Optimized system transfer function. 


Areay containing the phase cf the 
Simulation of the reed ; 
Optimized system transfer function. 


sp mnaag period if a zero crder 
hold is present in the system. 


Array ccntaining the discrete 
freguencies specified by the 
designer. 


Input Speed the maxinun 
frequency value of the range cf 
frequencies being ccnsideréd. 


Input specifying the mininun 
frequency value of the range cf 
frequencies being ccnsidered. 


Output giving the sampling 
frequency if a zero crder hold is 
present in the system. 


Array in which values to be varied 
by the minimizaticn routine are 
stored. At the start this aco 
contains the initial ques of the 
coefficient values and at . 
completion the ered contains the 
coefficient values that yield the 
Minimum cost functicn. 
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